076fecdb61815723d9bc72963e1a83d6334f5d3d
[u/mrichter/AliRoot.git] / test / vmctest / production / ConfigG3.C
1 // $Id$
2 //
3 // Configuration for the Geant3 simulations as Reference to Geant4 simulations
4 //
5
6 #if !defined(__CINT__) || defined(__MAKECINT__)
7 #include <Riostream.h>
8 #include <TRandom.h>
9 #include <TDatime.h>
10 #include <TSystem.h>
11 #include <TVirtualMC.h>
12 #include <TGeant3TGeo.h>
13 #include "STEER/AliRunLoader.h"
14 #include "STEER/AliRun.h"
15 #include "STEER/AliConfig.h"
16 #include "PYTHIA6/AliDecayerPythia.h"
17 #include "PYTHIA6/AliGenPythia.h"
18 #include "TDPMjet/AliGenDPMjet.h"
19 #include "STEER/AliMagFCheb.h"
20 #include "STRUCT/AliBODY.h"
21 #include "STRUCT/AliMAG.h"
22 #include "STRUCT/AliABSOv3.h"
23 #include "STRUCT/AliDIPOv3.h"
24 #include "STRUCT/AliHALLv3.h"
25 #include "STRUCT/AliFRAMEv2.h"
26 #include "STRUCT/AliSHILv3.h"
27 #include "STRUCT/AliPIPEv3.h"
28 #include "ITS/AliITSv11Hybrid.h"
29 #include "TPC/AliTPCv2.h"
30 #include "TOF/AliTOFv6T0.h"
31 #include "HMPID/AliHMPIDv3.h"
32 #include "ZDC/AliZDCv3.h"
33 #include "TRD/AliTRDv1.h"
34 #include "TRD/AliTRDgeometry.h"
35 #include "FMD/AliFMDv1.h"
36 #include "MUON/AliMUONv1.h"
37 #include "PHOS/AliPHOSv1.h"
38 #include "PHOS/AliPHOSSimParam.h"
39 #include "PMD/AliPMDv1.h"
40 #include "T0/AliT0v1.h"
41 #include "EMCAL/AliEMCALv2.h"
42 #include "ACORDE/AliACORDEv1.h"
43 #include "VZERO/AliVZEROv7.h"
44 #endif
45
46
47 enum PDC06Proc_t 
48   {
49     kPythia6, kPythia6D6T, kPythia6ATLAS, kPythia6ATLAS_Flat, kPythiaPerugia0, kPhojet, kRunMax
50   };
51
52 const char * pprRunName[] = {
53   "kPythia6", "kPythia6D6T", "kPythia6ATLAS", "kPythia6ATLAS_Flat", "kPythiaPerugia0", "kPhojet" 
54 };
55
56 enum Mag_t
57   {
58     kNoField, k5kG, kFieldMax
59   };
60
61 const char * pprField[] = {
62   "kNoField", "k5kG"
63 };
64
65
66 enum PprTrigConf_t
67   {
68     kDefaultPPTrig, kDefaultPbPbTrig
69   };
70
71 const char * pprTrigConfName[] = {
72   "p-p","Pb-Pb"
73 };
74
75 //--- Functions ---
76 class AliGenPythia;
77 AliGenerator *MbPythia();
78 AliGenerator *MbPythiaTuneD6T();
79 AliGenerator *MbPhojet();
80 void ProcessEnvironmentVars();
81
82 // Geterator, field, beam energy
83 static PDC06Proc_t   proc     = kPhojet;
84 static Mag_t         mag      = k5kG;
85 static Float_t       energy   = 10000; // energy in CMS
86 static Int_t         runNumber = 0;
87 static PprTrigConf_t strig = kDefaultPPTrig; // default pp trigger configuration
88 //========================//
89 // Set Random Number seed //
90 //========================//
91 TDatime dt;
92 static UInt_t seed    = dt.Get();
93
94 // Comment line
95 static TString comment;
96
97 void Config()
98 {
99     
100
101   // Get settings from environment variables
102   ProcessEnvironmentVars();
103
104   gRandom->SetSeed(seed);
105   cerr<<"Seed for random number generation= "<<seed<<endl; 
106
107   // Libraries required by geant321
108 #if defined(__CINT__)
109   gSystem->Load("liblhapdf");      // Parton density functions
110   gSystem->Load("libEGPythia6");   // TGenerator interface
111   if (proc == kPythia6 || proc == kPhojet) {
112     gSystem->Load("libpythia6");        // Pythia 6.2
113   } else {
114     gSystem->Load("libpythia6.4.21");   // Pythia 6.4
115   }
116   gSystem->Load("libAliPythia6");  // ALICE specific implementations
117   gSystem->Load("libgeant321");
118 #endif
119
120   new TGeant3TGeo("C++ Interface to Geant3");
121
122   //=======================================================================
123   //  Create the output file
124
125    
126   AliRunLoader* rl=0x0;
127
128   cout<<"Config.C: Creating Run Loader ..."<<endl;
129   rl = AliRunLoader::Open("galice.root",
130                           AliConfig::GetDefaultEventFolderName(),
131                           "recreate");
132   if (rl == 0x0)
133     {
134       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
135       return;
136     }
137   rl->SetCompressionLevel(2);
138   rl->SetNumberOfEventsPerFile(1000);
139   gAlice->SetRunLoader(rl);
140   // gAlice->SetGeometryFromFile("geometry.root");
141   // gAlice->SetGeometryFromCDB();
142   
143   // Set the trigger configuration: proton-proton
144
145   AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
146   cout <<"Trigger configuration is set to  "<<pprTrigConfName[strig]<<endl;
147   //
148   //=======================================================================
149   // ************* STEERING parameters FOR ALICE SIMULATION **************
150   // --- Specify event type to be tracked through the ALICE setup
151   // --- All positions are in cm, angles in degrees, and P and E in GeV
152
153
154   gMC->SetProcess("DCAY",1);
155   gMC->SetProcess("PAIR",1);
156   gMC->SetProcess("COMP",1);
157   gMC->SetProcess("PHOT",1);
158   gMC->SetProcess("PFIS",0);
159   gMC->SetProcess("DRAY",0);
160   gMC->SetProcess("ANNI",1);
161   gMC->SetProcess("BREM",1);
162   gMC->SetProcess("MUNU",1);
163   gMC->SetProcess("CKOV",1);
164   gMC->SetProcess("HADR",1);
165   gMC->SetProcess("LOSS",2);
166   gMC->SetProcess("MULS",1);
167   gMC->SetProcess("RAYL",1);
168
169   Float_t cut = 1.e-3;        // 1MeV cut by default
170   Float_t tofmax = 1.e10;
171
172   gMC->SetCut("CUTGAM", cut);
173   gMC->SetCut("CUTELE", cut);
174   gMC->SetCut("CUTNEU", cut);
175   gMC->SetCut("CUTHAD", cut);
176   gMC->SetCut("CUTMUO", cut);
177   gMC->SetCut("BCUTE",  cut); 
178   gMC->SetCut("BCUTM",  cut); 
179   gMC->SetCut("DCUTE",  cut); 
180   gMC->SetCut("DCUTM",  cut); 
181   gMC->SetCut("PPCUTM", cut);
182   gMC->SetCut("TOFMAX", tofmax); 
183
184
185
186
187   //======================//
188   // Set External decayer //
189   //======================//
190   TVirtualMCDecayer* decayer = new AliDecayerPythia();
191   decayer->SetForceDecay(kAll);
192   decayer->Init();
193   gMC->SetExternalDecayer(decayer);
194
195   //=========================//
196   // Generator Configuration //
197   //=========================//
198   AliGenerator* gener = 0x0;
199   
200   if (proc == kPythia6) {
201     gener = MbPythia();
202   } else if (proc == kPythia6D6T) {
203     gener = MbPythiaTuneD6T();
204   } else if (proc == kPythia6ATLAS) {
205     gener = MbPythiaTuneATLAS();
206   } else if (proc == kPythiaPerugia0) {
207     gener = MbPythiaTunePerugia0();
208   } else if (proc == kPythia6ATLAS_Flat) {
209     gener = MbPythiaTuneATLAS_Flat();
210   } else if (proc == kPhojet) {
211     gener = MbPhojet();
212   }
213   
214   
215   //
216   //
217   // Size of the interaction diamond
218   // Longitudinal
219   Float_t sigmaz  = 5.4 / TMath::Sqrt(2.); // [cm]
220   if (energy == 900)
221     //sigmaz  = 10.5 / TMath::Sqrt(2.); // [cm]
222     //sigmaz = 3.7;
223     if (energy == 7000)
224       sigmaz  = 6.3 / TMath::Sqrt(2.); // [cm]
225   
226   //
227   // Transverse
228
229   // beta*
230   Float_t betast                  = 10.0;  // beta* [m]
231   if (runNumber >= 117048) betast =  2.0;
232   if (runNumber >  122375) betast =  3.5;  // starting with fill 1179
233   //    
234   //
235   Float_t eps     = 5.0e-6;            // emittance [m]
236   Float_t gamma   = energy / 2.0 / 0.938272;  // relativistic gamma [1]
237   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.;  // [cm]
238   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
239     
240   gener->SetSigma(sigmaxy, sigmaxy, sigmaz);      // Sigma in (X,Y,Z) (cm) on IP position
241   gener->SetVertexSmear(kPerEvent);
242   gener->Init();
243
244   printf("\n \n Comment: %s \n \n", comment.Data());
245
246   rl->CdGAFile();
247   
248   Int_t iABSO  = 1;
249   Int_t iACORDE= 0;
250   Int_t iDIPO  = 1;
251   Int_t iEMCAL = 1;
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   = 1;
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
271   //=================== Alice BODY parameters =============================
272   AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
273
274
275   if (iMAG)
276     {
277       //=================== MAG parameters ============================
278       // --- Start with Magnet since detector layouts may be depending ---
279       // --- on the selected Magnet dimensions ---
280       AliMAG *MAG = new AliMAG("MAG", "Magnet");
281     }
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
305   if (iFRAME)
306     {
307       //=================== FRAME parameters ============================
308
309       AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
310       FRAME->SetHoles(1);
311     }
312
313   if (iSHIL)
314     {
315       //=================== SHIL parameters ============================
316
317       AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
318     }
319
320
321   if (iPIPE)
322     {
323       //=================== PIPE parameters ============================
324
325       AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
326     }
327  
328   if (iITS)
329     {
330       //=================== ITS parameters ============================
331
332       //AliITS *ITS  = new AliITSv11Hybrid("ITS","ITS v11Hybrid");
333       AliITS *ITS  = new AliITSv11("ITS","ITS v11");
334     }
335
336   if (iTPC)
337     {
338       //============================ TPC parameters =====================
339
340       AliTPC *TPC = new AliTPCv2("TPC", "Default");
341     }
342
343
344   if (iTOF) {
345     //=================== TOF parameters ============================
346
347     AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
348   }
349
350
351   if (iHMPID)
352     {
353       //=================== HMPID parameters ===========================
354
355       AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
356
357     }
358
359
360   if (iZDC)
361     {
362       //=================== ZDC parameters ============================
363         
364       AliZDC *ZDC = new AliZDCv3("ZDC", "normal ZDC");
365       //Collimators aperture
366       ZDC->SetVCollSideCAperture(0.85);
367       ZDC->SetVCollSideCCentre(0.);
368       ZDC->SetVCollSideAAperture(0.75);
369       ZDC->SetVCollSideACentre(0.);
370       //Detector position
371       ZDC->SetYZNC(1.6);
372       ZDC->SetYZNA(1.6);
373       ZDC->SetYZPC(1.6);
374       ZDC->SetYZPA(1.6);
375     }
376
377   if (iTRD)
378     {
379       //=================== TRD parameters ============================
380
381       AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
382       AliTRDgeometry *geoTRD = TRD->GetGeometry();
383       // Partial geometry: modules at 0,1,7,8,9,16,17
384       // starting at 3h in positive direction
385       geoTRD->SetSMstatus(2,0);
386       geoTRD->SetSMstatus(3,0);
387       geoTRD->SetSMstatus(4,0);
388       geoTRD->SetSMstatus(5,0);
389       geoTRD->SetSMstatus(6,0);
390       geoTRD->SetSMstatus(11,0);
391       geoTRD->SetSMstatus(12,0);
392       geoTRD->SetSMstatus(13,0);
393       geoTRD->SetSMstatus(14,0);
394       geoTRD->SetSMstatus(15,0);
395       geoTRD->SetSMstatus(16,0);
396     }
397
398   if (iFMD)
399     {
400       //=================== FMD parameters ============================
401
402       AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
403     }
404
405   if (iMUON)
406     {
407       //=================== MUON parameters ===========================
408       // New MUONv1 version (geometry defined via builders)
409       AliMUON *MUON = new AliMUONv1("MUON", "default");
410       // activate trigger efficiency by cells
411       MUON->SetTriggerEffCells(1);
412     }
413
414   if (iPHOS)
415     {
416       //=================== PHOS parameters ===========================
417
418       AliPHOS *PHOS = new AliPHOSv1("PHOS", "noCPV_Modules123");
419
420     }
421
422
423   if (iPMD)
424     {
425       //=================== PMD parameters ============================
426
427       AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
428     }
429
430   if (iT0)
431     {
432       //=================== T0 parameters ============================
433       AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
434     }
435
436   if (iEMCAL)
437     {
438       //=================== EMCAL parameters ============================
439
440       AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_FIRSTYEARV1");
441     }
442
443   if (iACORDE)
444     {
445       //=================== ACORDE parameters ============================
446
447       AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
448     }
449
450   if (iVZERO)
451     {
452       //=================== ACORDE parameters ============================
453
454       AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
455     }
456 }
457 //
458 //           PYTHIA
459 //
460
461 AliGenerator* MbPythia()
462 {
463   comment = comment.Append(" pp: Pythia low-pt");
464   //
465   //    Pythia
466   AliGenPythia* pythia = new AliGenPythia(-1); 
467   pythia->SetMomentumRange(0, 999999.);
468   pythia->SetThetaRange(0., 180.);
469   pythia->SetYRange(-12.,12.);
470   pythia->SetPtRange(0,1000.);
471   pythia->SetProcess(kPyMb);
472   pythia->SetEnergyCMS(energy);
473       
474   return pythia;
475 }
476
477 AliGenerator* MbPythiaTuneD6T()
478 {
479   comment = comment.Append(" pp: Pythia low-pt");
480   //
481   //    Pythia
482   AliGenPythia* pythia = new AliGenPythia(-1); 
483   pythia->SetMomentumRange(0, 999999.);
484   pythia->SetThetaRange(0., 180.);
485   pythia->SetYRange(-12.,12.);
486   pythia->SetPtRange(0,1000.);
487   pythia->SetProcess(kPyMb);
488   pythia->SetEnergyCMS(energy);
489   //    Tune
490   //    109     D6T : Rick Field's CDF Tune D6T (NB: needs CTEQ6L pdfs externally)
491   pythia->SetTune(109); // F I X 
492   pythia->SetStrucFunc(kCTEQ6l);
493   //
494   return pythia;
495 }
496
497 AliGenerator* MbPythiaTunePerugia0()
498 {
499   comment = comment.Append(" pp: Pythia low-pt (Perugia0)");
500   //
501   //    Pythia
502   AliGenPythia* pythia = new AliGenPythia(-1); 
503   pythia->SetMomentumRange(0, 999999.);
504   pythia->SetThetaRange(0., 180.);
505   pythia->SetYRange(-12.,12.);
506   pythia->SetPtRange(0,1000.);
507   pythia->SetProcess(kPyMb);
508   pythia->SetEnergyCMS(energy);
509   //    Tune
510   //    320     Perugia 0
511   pythia->SetTune(320); 
512   pythia->UseNewMultipleInteractionsScenario();
513   pythia->SetCrossingAngle(0,0.000280);
514
515   //
516   return pythia;
517 }
518
519
520 AliGenerator* MbPythiaTuneATLAS()
521 {
522   comment = comment.Append(" pp: Pythia low-pt");
523   //
524   //    Pythia
525   AliGenPythia* pythia = new AliGenPythia(-1); 
526   pythia->SetMomentumRange(0, 999999.);
527   pythia->SetThetaRange(0., 180.);
528   pythia->SetYRange(-12.,12.);
529   pythia->SetPtRange(0,1000.);
530   pythia->SetProcess(kPyMb);
531   pythia->SetEnergyCMS(energy);
532   //    Tune
533   //    C   306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune (needs CTEQ6L externally)
534   pythia->SetTune(306);
535   pythia->SetStrucFunc(kCTEQ6l);
536   //
537   return pythia;
538 }
539
540 AliGenerator* MbPythiaTuneATLAS_Flat()
541 {
542   AliGenPythia* pythia = MbPythiaTuneATLAS();
543       
544   comment = comment.Append("; flat multiplicity distribution");
545       
546   // set high multiplicity trigger
547   // this weight achieves a flat multiplicity distribution
548   Double_t cont[] =
549     {0, 
550      0.234836, 0.103484, 0.00984802, 0.0199906, 0.0260018, 0.0208481, 0.0101477, 0.00146998, -0.00681513, -0.0114928,
551      -0.0113352, -0.0102012, -0.00895238, -0.00534961, -0.00261648, -0.000819048, 0.00130816, 0.00177978, 0.00373838, 0.00566255,
552      0.00628156, 0.00687458, 0.00668017, 0.00702917, 0.00810848, 0.00876167, 0.00832783, 0.00848518, 0.0107709, 0.0106849,
553      0.00933702, 0.00912525, 0.0106553, 0.0102785, 0.0101756, 0.010962, 0.00957103, 0.00970448, 0.0117133, 0.012271,
554      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
555      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
556      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
557      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
558      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
559      0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113, 0.0113,
560      0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
561      0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
562      0};
563
564   Double_t err[] =
565     {0, 
566      0.00168216, 0.000743548, 0.00103125, 0.00108605, 0.00117101, 0.00124577, 0.00129119, 0.00128341, 0.00121421, 0.00112431,
567      0.00100588, 0.000895078, 0.000790314, 0.000711673, 0.000634547, 0.000589133, 0.000542763, 0.000509552, 0.000487375, 0.000468906, 
568      0.000460196, 0.000455806, 0.00044843, 0.000449317, 0.00045007, 0.000458016, 0.000461167, 0.000474742, 0.00050161, 0.00051637, 
569      0.000542336, 0.000558854, 0.000599169, 0.000611982, 0.000663855, 0.000696322, 0.000722825, 0.000771686, 0.000838023, 0.000908317, 
570      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
571      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
572      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
573      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
574      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
575      0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003,
576      0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
577      0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
578      0};
579
580   TH1F *weight = new TH1F("newweight","newweight",120,-0.5,119.5);
581
582   weight->SetContent(cont);
583   weight->SetError(err);
584         
585   Int_t limit = weight->GetRandom();
586   pythia->SetTriggerChargedMultiplicity(limit, 1.4);
587       
588   comment = comment.Append(Form("; multiplicity threshold set to %d in |eta| < 1.4", limit));
589
590   return pythia;
591 }
592
593 AliGenerator* MbPhojet()
594 {
595   comment = comment.Append(" pp: Pythia low-pt");
596   //
597   //    DPMJET
598 #if defined(__CINT__)
599   gSystem->Load("libdpmjet");      // Parton density functions
600   gSystem->Load("libTDPMjet");      // Parton density functions
601 #endif
602   AliGenDPMjet* dpmjet = new AliGenDPMjet(-1); 
603   dpmjet->SetMomentumRange(0, 999999.);
604   dpmjet->SetThetaRange(0., 180.);
605   dpmjet->SetYRange(-12.,12.);
606   dpmjet->SetPtRange(0,1000.);
607   dpmjet->SetProcess(kDpmMb);
608   dpmjet->SetEnergyCMS(energy);
609   dpmjet->SetCrossingAngle(0,0.000280);
610   return dpmjet;
611 }
612
613 void ProcessEnvironmentVars()
614 {
615   // Run type
616   if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
617     for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
618       if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0) {
619         proc = (PDC06Proc_t)iRun;
620         cout<<"Run type set to "<<pprRunName[iRun]<<endl;
621       }
622     }
623   }
624
625   // Field
626   if (gSystem->Getenv("CONFIG_FIELD")) {
627     for (Int_t iField = 0; iField < kFieldMax; iField++) {
628       if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField])==0) {
629         mag = (Mag_t)iField;
630         cout<<"Field set to "<<pprField[iField]<<endl;
631       }
632     }
633   }
634
635   // Energy
636   if (gSystem->Getenv("CONFIG_ENERGY")) {
637     energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
638     cout<<"Energy set to "<<energy<<" GeV"<<endl;
639   }
640
641   // Random Number seed
642   if (gSystem->Getenv("CONFIG_SEED")) {
643     seed = atoi(gSystem->Getenv("CONFIG_SEED"));
644   }
645
646   // Run number
647   if (gSystem->Getenv("DC_RUN")) {
648     runNumber = atoi(gSystem->Getenv("DC_RUN"));
649   }
650 }