f0eb2f54bc9eae6fbf9644048fced0018cca44c7
[u/mrichter/AliRoot.git] / PWGGA / PHOSTasks / PHOS_embedding / 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/AliITSv11Hybrid.h"
34 #include "TPC/AliTPCv2.h"
35 #include "TOF/AliTOFv6T0.h"
36 #include "HMPID/AliHMPIDv3.h"
37 #include "ZDC/AliZDCv3.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
52 enum PDC06Proc_t 
53 {
54   kPythia6, kPythia6D6T, kPythia6ATLAS, kPythia6ATLAS_Flat, kPythiaPerugia0, kPhojet, kRunMax
55 };
56
57 const char * pprRunName[] = {
58   "kPythia6", "kPythia6D6T", "kPythia6ATLAS", "kPythia6ATLAS_Flat", "kPythiaPerugia0", "kPhojet" 
59 };
60
61 enum Mag_t
62 {
63   kNoField, k5kG, kFieldMax
64 };
65
66 const char * pprField[] = {
67   "kNoField", "k5kG"
68 };
69
70 //--- Functions ---
71 class AliGenPythia;
72 AliGenerator *MbPythia();
73 AliGenerator *MbPythiaTuneD6T();
74 AliGenerator *MbPhojet();
75 void ProcessEnvironmentVars();
76
77 // Geterator, field, beam energy
78 static PDC06Proc_t   proc     = kPhojet;
79 static Mag_t         mag      = k5kG;
80 static Float_t       energy   = 10000; // energy in CMS
81 static Int_t         runNumber = 0;
82
83 //========================//
84 // Set Random Number seed //
85 //========================//
86 TDatime dt;
87 static UInt_t seed    = dt.Get();
88
89 // Comment line
90 static TString comment;
91
92 void Config()
93 {
94     
95
96   // Get settings from environment variables
97   ProcessEnvironmentVars();
98
99   gRandom->SetSeed(seed);
100   cerr<<"Seed for random number generation= "<<seed<<endl; 
101
102   // Libraries required by geant321
103 #if defined(__CINT__)
104   gSystem->Load("liblhapdf");      // Parton density functions
105   gSystem->Load("libEGPythia6");   // TGenerator interface
106   if (proc == kPythia6 || proc == kPhojet) {
107     gSystem->Load("libpythia6");        // Pythia 6.2
108   } else {
109     gSystem->Load("libpythia6.4.21");   // Pythia 6.4
110   }
111   gSystem->Load("libAliPythia6");  // ALICE specific implementations
112   gSystem->Load("libgeant321");
113 #endif
114
115   new TGeant3TGeo("C++ Interface to Geant3");
116
117   //=======================================================================
118   //  Create the output file
119
120    
121   AliRunLoader* rl=0x0;
122
123   cout<<"Config.C: Creating Run Loader ..."<<endl;
124   rl = AliRunLoader::Open("galice.root",
125                           AliConfig::GetDefaultEventFolderName(),
126                           "recreate");
127   if (rl == 0x0)
128     {
129       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
130       return;
131     }
132   rl->SetCompressionLevel(2);
133   rl->SetNumberOfEventsPerFile(3000);
134   gAlice->SetRunLoader(rl);
135   // gAlice->SetGeometryFromFile("geometry.root");
136   // gAlice->SetGeometryFromCDB();
137   
138   // Set the trigger configuration: proton-proton
139  // gAlice->SetTriggerDescriptor("p-p");
140
141   //
142   //=======================================================================
143   // ************* STEERING parameters FOR ALICE SIMULATION **************
144   // --- Specify event type to be tracked through the ALICE setup
145   // --- All positions are in cm, angles in degrees, and P and E in GeV
146
147
148     gMC->SetProcess("DCAY",1);
149     gMC->SetProcess("PAIR",1);
150     gMC->SetProcess("COMP",1);
151     gMC->SetProcess("PHOT",1);
152     gMC->SetProcess("PFIS",0);
153     gMC->SetProcess("DRAY",0);
154     gMC->SetProcess("ANNI",1);
155     gMC->SetProcess("BREM",1);
156     gMC->SetProcess("MUNU",1);
157     gMC->SetProcess("CKOV",1);
158     gMC->SetProcess("HADR",1);
159     gMC->SetProcess("LOSS",2);
160     gMC->SetProcess("MULS",1);
161     gMC->SetProcess("RAYL",1);
162
163     Float_t cut = 1.e-3;        // 1MeV cut by default
164     Float_t tofmax = 1.e10;
165
166     gMC->SetCut("CUTGAM", cut);
167     gMC->SetCut("CUTELE", cut);
168     gMC->SetCut("CUTNEU", cut);
169     gMC->SetCut("CUTHAD", cut);
170     gMC->SetCut("CUTMUO", cut);
171     gMC->SetCut("BCUTE",  cut); 
172     gMC->SetCut("BCUTM",  cut); 
173     gMC->SetCut("DCUTE",  cut); 
174     gMC->SetCut("DCUTM",  cut); 
175     gMC->SetCut("PPCUTM", cut);
176     gMC->SetCut("TOFMAX", tofmax); 
177
178
179
180
181   //======================//
182   // Set External decayer //
183   //======================//
184   TVirtualMCDecayer* decayer = new AliDecayerPythia();
185   decayer->SetForceDecay(kAll);
186   decayer->Init();
187   gMC->SetExternalDecayer(decayer);
188
189   //=========================//
190   // Generator Configuration //
191   //=========================//
192 /*
193   AliGenerator* gener = 0x0;
194   
195   if (proc == kPythia6) {
196       gener = MbPythia();
197   } else if (proc == kPythia6D6T) {
198       gener = MbPythiaTuneD6T();
199   } else if (proc == kPythia6ATLAS) {
200       gener = MbPythiaTuneATLAS();
201   } else if (proc == kPythiaPerugia0) {
202       gener = MbPythiaTunePerugia0();
203   } else if (proc == kPythia6ATLAS_Flat) {
204       gener = MbPythiaTuneATLAS_Flat();
205   } else if (proc == kPhojet) {
206       gener = MbPhojet();
207   }
208 */  
209 /*
210    AliGenCocktail *gener  = new AliGenCocktail();
211
212    AliGenPHOSlib * lib = new AliGenPHOSlib() ;
213    //4pi0
214    AliGenParam *genPHOS = new AliGenParam(4,AliGenPHOSlib::kPion,lib->GetPt(AliGenPHOSlib::kPi0),lib->GetY(AliGenPHOSlib::kPi0Flat),lib->GetIp(AliGenPHOSlib::kPi0)) ;
215    genPHOS->SetPhiRange(250.,330.) ;
216    genPHOS->SetYRange(-0.15.,0.15) ;
217    genPHOS->SetPtRange(0.5,30.) ;
218    gener->AddGenerator(genPHOS,"PHOS",1.) ;
219 */
220
221 printf("Creating FLAT generator \n") ;
222
223    gener = new AliGenBox(2) ;
224    gener->SetPhiRange(250.,330.) ;
225    gener->SetYRange(-0.15.,0.15) ;
226    gener->SetEtaRange(-0.15,0.15) ;
227    gener->SetPtRange(0.5,30.) ;
228    gener->SetPart(111) ;
229
230   //
231   //
232   // Size of the interaction diamond
233   // Longitudinal
234   Float_t sigmaz  = 5.4 / TMath::Sqrt(2.); // [cm]
235   if (energy == 900)
236     //sigmaz  = 10.5 / TMath::Sqrt(2.); // [cm]
237     //sigmaz = 3.7;
238   if (energy == 7000)
239     sigmaz  = 6.3 / TMath::Sqrt(2.); // [cm]
240   
241   //
242   // Transverse
243   // beta*
244   Float_t betast  = 10.;                 // beta* [m]
245   if (runNumber >= 117048) betast = 2.;
246   printf("beta* for run# %8d is %13.3f", runNumber, betast);
247   //
248   Float_t eps     = 3.75e-6;            // emittance [m]
249   Float_t gamma   = energy / 2.0 / 0.938272;  // relativistic gamma [1]
250   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.;  // [cm]
251   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
252     
253   gener->SetSigma(sigmaxy, sigmaxy, sigmaz);      // Sigma in (X,Y,Z) (cm) on IP position
254   gener->SetVertexSmear(kPerEvent);
255   gener->Init();
256
257   printf("\n \n Comment: %s \n \n", comment.Data());
258
259   rl->CdGAFile();
260   
261   Int_t iABSO  = 1;
262   Int_t iACORDE= 0;
263   Int_t iDIPO  = 1;
264   Int_t iEMCAL = 1;
265   Int_t iFMD   = 1;
266   Int_t iFRAME = 1;
267   Int_t iHALL  = 1;
268   Int_t iITS   = 1;
269   Int_t iMAG   = 1;
270   Int_t iMUON  = 1;
271   Int_t iPHOS  = 1;
272   Int_t iPIPE  = 1;
273   Int_t iPMD   = 1;
274   Int_t iHMPID = 1;
275   Int_t iSHIL  = 1;
276   Int_t iT0    = 1;
277   Int_t iTOF   = 1;
278   Int_t iTPC   = 1;
279   Int_t iTRD   = 1;
280   Int_t iVZERO = 1;
281   Int_t iZDC   = 1;
282   
283
284     //=================== Alice BODY parameters =============================
285     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
286
287
288     if (iMAG)
289     {
290         //=================== MAG parameters ============================
291         // --- Start with Magnet since detector layouts may be depending ---
292         // --- on the selected Magnet dimensions ---
293         AliMAG *MAG = new AliMAG("MAG", "Magnet");
294     }
295
296
297     if (iABSO)
298     {
299         //=================== ABSO parameters ============================
300         AliABSO *ABSO = new AliABSOv3("ABSO", "Muon Absorber");
301     }
302
303     if (iDIPO)
304     {
305         //=================== DIPO parameters ============================
306
307         AliDIPO *DIPO = new AliDIPOv3("DIPO", "Dipole version 3");
308     }
309
310     if (iHALL)
311     {
312         //=================== HALL parameters ============================
313
314         AliHALL *HALL = new AliHALLv3("HALL", "Alice Hall");
315     }
316
317
318     if (iFRAME)
319     {
320         //=================== FRAME parameters ============================
321
322         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
323         FRAME->SetHoles(1);
324     }
325
326     if (iSHIL)
327     {
328         //=================== SHIL parameters ============================
329
330         AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
331     }
332
333
334     if (iPIPE)
335     {
336         //=================== PIPE parameters ============================
337
338         AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
339     }
340  
341     if (iITS)
342     {
343         //=================== ITS parameters ============================
344
345         AliITS *ITS  = new AliITSv11Hybrid("ITS","ITS v11Hybrid");
346     }
347
348     if (iTPC)
349     {
350       //============================ TPC parameters =====================
351
352         AliTPC *TPC = new AliTPCv2("TPC", "Default");
353     }
354
355
356     if (iTOF) {
357         //=================== TOF parameters ============================
358
359         AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
360     }
361
362
363     if (iHMPID)
364     {
365         //=================== HMPID parameters ===========================
366
367         AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
368
369     }
370
371
372     if (iZDC)
373     {
374         //=================== ZDC parameters ============================
375
376         AliZDC *ZDC = new AliZDCv3("ZDC", "normal ZDC");
377     }
378
379     if (iTRD)
380     {
381         //=================== TRD parameters ============================
382
383         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
384         AliTRDgeometry *geoTRD = TRD->GetGeometry();
385         // Partial geometry: modules at 0,1,7,8,9,16,17
386         // starting at 3h in positive direction
387         geoTRD->SetSMstatus(2,0);
388         geoTRD->SetSMstatus(3,0);
389         geoTRD->SetSMstatus(4,0);
390         geoTRD->SetSMstatus(5,0);
391         geoTRD->SetSMstatus(6,0);
392         geoTRD->SetSMstatus(11,0);
393         geoTRD->SetSMstatus(12,0);
394         geoTRD->SetSMstatus(13,0);
395         geoTRD->SetSMstatus(14,0);
396         geoTRD->SetSMstatus(15,0);
397         geoTRD->SetSMstatus(16,0);
398     }
399
400     if (iFMD)
401     {
402         //=================== FMD parameters ============================
403
404         AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
405    }
406
407     if (iMUON)
408     {
409         //=================== MUON parameters ===========================
410         // New MUONv1 version (geometry defined via builders)
411
412         AliMUON *MUON = new AliMUONv1("MUON", "default");
413         
414         // activate trigger efficiency by cells
415
416         MUON->SetTriggerEffCells(1);     
417
418      }
419
420     if (iPHOS)
421     {
422         //=================== PHOS parameters ===========================
423
424      AliPHOS *PHOS = new AliPHOSv1("PHOS", "noCPV_Modules123");
425
426     }
427
428
429     if (iPMD)
430     {
431         //=================== PMD parameters ============================
432
433         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
434     }
435
436     if (iT0)
437     {
438         //=================== T0 parameters ============================
439         AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
440     }
441
442     if (iEMCAL)
443     {
444         //=================== EMCAL parameters ============================
445
446         AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_FIRSTYEAR");
447     }
448
449      if (iACORDE)
450     {
451         //=================== ACORDE parameters ============================
452
453         AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
454     }
455
456      if (iVZERO)
457     {
458         //=================== ACORDE parameters ============================
459
460         AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
461     }
462 }
463 //
464 //           PYTHIA
465 //
466
467 AliGenerator* MbPythia()
468 {
469       comment = comment.Append(" pp: Pythia low-pt");
470 //
471 //    Pythia
472       AliGenPythia* pythia = new AliGenPythia(-1); 
473       pythia->SetMomentumRange(0, 999999.);
474       pythia->SetThetaRange(0., 180.);
475       pythia->SetYRange(-12.,12.);
476       pythia->SetPtRange(0,1000.);
477       pythia->SetProcess(kPyMb);
478       pythia->SetEnergyCMS(energy);
479       
480       return pythia;
481 }
482
483 AliGenerator* MbPythiaTuneD6T()
484 {
485       comment = comment.Append(" pp: Pythia low-pt");
486 //
487 //    Pythia
488       AliGenPythia* pythia = new AliGenPythia(-1); 
489       pythia->SetMomentumRange(0, 999999.);
490       pythia->SetThetaRange(0., 180.);
491       pythia->SetYRange(-12.,12.);
492       pythia->SetPtRange(0,1000.);
493       pythia->SetProcess(kPyMb);
494       pythia->SetEnergyCMS(energy);
495 //    Tune
496 //    109     D6T : Rick Field's CDF Tune D6T (NB: needs CTEQ6L pdfs externally)
497       pythia->SetTune(109); // F I X 
498       pythia->SetStrucFunc(kCTEQ6l);
499 //
500       return pythia;
501 }
502
503 AliGenerator* MbPythiaTunePerugia0()
504 {
505       comment = comment.Append(" pp: Pythia low-pt (Perugia0)");
506 //
507 //    Pythia
508       AliGenPythia* pythia = new AliGenPythia(-1); 
509       pythia->SetMomentumRange(0, 999999.);
510       pythia->SetThetaRange(0., 180.);
511       pythia->SetYRange(-12.,12.);
512       pythia->SetPtRange(0,1000.);
513       pythia->SetProcess(kPyMb);
514       pythia->SetEnergyCMS(energy);
515 //    Tune
516 //    320     Perugia 0
517       pythia->SetTune(320); 
518       pythia->UseNewMultipleInteractionsScenario();
519 //
520       return pythia;
521 }
522
523
524 AliGenerator* MbPythiaTuneATLAS()
525 {
526       comment = comment.Append(" pp: Pythia low-pt");
527 //
528 //    Pythia
529       AliGenPythia* pythia = new AliGenPythia(-1); 
530       pythia->SetMomentumRange(0, 999999.);
531       pythia->SetThetaRange(0., 180.);
532       pythia->SetYRange(-12.,12.);
533       pythia->SetPtRange(0,1000.);
534       pythia->SetProcess(kPyMb);
535       pythia->SetEnergyCMS(energy);
536 //    Tune
537 //    C   306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune (needs CTEQ6L externally)
538       pythia->SetTune(306);
539       pythia->SetStrucFunc(kCTEQ6l);
540 //
541       return pythia;
542 }
543
544 AliGenerator* MbPythiaTuneATLAS_Flat()
545 {
546       AliGenPythia* pythia = MbPythiaTuneATLAS();
547       
548       comment = comment.Append("; flat multiplicity distribution");
549       
550       // set high multiplicity trigger
551       // this weight achieves a flat multiplicity distribution
552       TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
553       weight->SetBinContent(1,5.49443);
554       weight->SetBinContent(2,8.770816);
555       weight->SetBinContent(6,0.4568624);
556       weight->SetBinContent(7,0.2919915);
557       weight->SetBinContent(8,0.6674189);
558       weight->SetBinContent(9,0.364737);
559       weight->SetBinContent(10,0.8818444);
560       weight->SetBinContent(11,0.531885);
561       weight->SetBinContent(12,1.035197);
562       weight->SetBinContent(13,0.9394057);
563       weight->SetBinContent(14,0.9643193);
564       weight->SetBinContent(15,0.94543);
565       weight->SetBinContent(16,0.9426507);
566       weight->SetBinContent(17,0.9423649);
567       weight->SetBinContent(18,0.789456);
568       weight->SetBinContent(19,1.149026);
569       weight->SetBinContent(20,1.100491);
570       weight->SetBinContent(21,0.6350525);
571       weight->SetBinContent(22,1.351941);
572       weight->SetBinContent(23,0.03233504);
573       weight->SetBinContent(24,0.9574557);
574       weight->SetBinContent(25,0.868133);
575       weight->SetBinContent(26,1.030998);
576       weight->SetBinContent(27,1.08897);
577       weight->SetBinContent(28,1.251382);
578       weight->SetBinContent(29,0.1391099);
579       weight->SetBinContent(30,1.192876);
580       weight->SetBinContent(31,0.448944);
581       weight->SetBinContent(32,1);
582       weight->SetBinContent(33,1);
583       weight->SetBinContent(34,1);
584       weight->SetBinContent(35,1);
585       weight->SetBinContent(36,0.9999997);
586       weight->SetBinContent(37,0.9999997);
587       weight->SetBinContent(38,0.9999996);
588       weight->SetBinContent(39,0.9999996);
589       weight->SetBinContent(40,0.9999995);
590       weight->SetBinContent(41,0.9999993);
591       weight->SetBinContent(42,1);
592       weight->SetBinContent(43,1);
593       weight->SetBinContent(44,1);
594       weight->SetBinContent(45,1);
595       weight->SetBinContent(46,1);
596       weight->SetBinContent(47,0.9999999);
597       weight->SetBinContent(48,0.9999998);
598       weight->SetBinContent(49,0.9999998);
599       weight->SetBinContent(50,0.9999999);
600       weight->SetBinContent(51,0.9999999);
601       weight->SetBinContent(52,0.9999999);
602       weight->SetBinContent(53,0.9999999);
603       weight->SetBinContent(54,0.9999998);
604       weight->SetBinContent(55,0.9999998);
605       weight->SetBinContent(56,0.9999998);
606       weight->SetBinContent(57,0.9999997);
607       weight->SetBinContent(58,0.9999996);
608       weight->SetBinContent(59,0.9999995);
609       weight->SetBinContent(60,1);
610       weight->SetBinContent(61,1);
611       weight->SetBinContent(62,1);
612       weight->SetBinContent(63,1);
613       weight->SetBinContent(64,1);
614       weight->SetBinContent(65,0.9999999);
615       weight->SetBinContent(66,0.9999998);
616       weight->SetBinContent(67,0.9999998);
617       weight->SetBinContent(68,0.9999999);
618       weight->SetBinContent(69,1);
619       weight->SetBinContent(70,1);
620       weight->SetBinContent(71,0.9999997);
621       weight->SetBinContent(72,0.9999995);
622       weight->SetBinContent(73,0.9999994);
623       weight->SetBinContent(74,1);
624       weight->SetBinContent(75,1);
625       weight->SetBinContent(76,1);
626       weight->SetBinContent(77,1);
627       weight->SetBinContent(78,0.9999999);
628       weight->SetBinContent(79,1);
629       weight->SetBinContent(80,1);
630       weight->SetEntries(526);
631         
632       Int_t limit = weight->GetRandom();
633       pythia->SetTriggerChargedMultiplicity(limit, 1.4);
634       
635       comment = comment.Append(Form("; multiplicity threshold set to %d in |eta| < 1.4", limit));
636
637       return pythia;
638 }
639
640 AliGenerator* MbPhojet()
641 {
642       comment = comment.Append(" pp: Pythia low-pt");
643 //
644 //    DPMJET
645 #if defined(__CINT__)
646   gSystem->Load("libdpmjet");      // Parton density functions
647   gSystem->Load("libTDPMjet");      // Parton density functions
648 #endif
649       AliGenDPMjet* dpmjet = new AliGenDPMjet(-1); 
650       dpmjet->SetMomentumRange(0, 999999.);
651       dpmjet->SetThetaRange(0., 180.);
652       dpmjet->SetYRange(-12.,12.);
653       dpmjet->SetPtRange(0,1000.);
654       dpmjet->SetProcess(kDpmMb);
655       dpmjet->SetEnergyCMS(energy);
656
657       return dpmjet;
658 }
659
660 void ProcessEnvironmentVars()
661 {
662     // Run type
663     if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
664       for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
665         if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0) {
666           proc = (PDC06Proc_t)iRun;
667           cout<<"Run type set to "<<pprRunName[iRun]<<endl;
668         }
669       }
670     }
671
672     // Field
673     if (gSystem->Getenv("CONFIG_FIELD")) {
674       for (Int_t iField = 0; iField < kFieldMax; iField++) {
675         if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField])==0) {
676           mag = (Mag_t)iField;
677           cout<<"Field set to "<<pprField[iField]<<endl;
678         }
679       }
680     }
681
682     // Energy
683     if (gSystem->Getenv("CONFIG_ENERGY")) {
684       energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
685       cout<<"Energy set to "<<energy<<" GeV"<<endl;
686     }
687
688     // Random Number seed
689     if (gSystem->Getenv("CONFIG_SEED")) {
690       seed = atoi(gSystem->Getenv("CONFIG_SEED"));
691     }
692
693     // Run number
694     if (gSystem->Getenv("DC_RUN")) {
695       runNumber = atoi(gSystem->Getenv("DC_RUN"));
696     }
697 }