Adding all jet cases. Returning to VZERO v2
[u/mrichter/AliRoot.git] / macros / ConfigPPR.C
1 // One can use the configuration macro in compiled mode by
2 // root [0] gSystem->Load("libgeant321");
3 // root [0] gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include\
4 //                   -I$ALICE_ROOT -I$ALICE/geant3/TGeant3");
5 // root [0] .x grun.C(1,"ConfigPPR.C++")
6
7 #if !defined(__CINT__) || defined(__MAKECINT__)
8 #include <Riostream.h>
9 #include <TRandom.h>
10 #include <TSystem.h>
11 #include <TVirtualMC.h>
12 #include <TGeant3.h>
13 #include "STEER/AliRunLoader.h"
14 #include "STEER/AliRun.h"
15 #include "STEER/AliConfig.h"
16 #include "STEER/AliGenerator.h"
17 #include "PYTHIA6/AliDecayerPythia.h"
18 #include "EVGEN/AliGenHIJINGpara.h"
19 #include "THijing/AliGenHijing.h"
20 #include "EVGEN/AliGenCocktail.h"
21 #include "EVGEN/AliGenSlowNucleons.h"
22 #include "EVGEN/AliSlowNucleonModelExp.h"
23 #include "EVGEN/AliGenParam.h"
24 #include "EVGEN/AliGenGSIlib.h"
25 #include "PYTHIA6/AliGenPythia.h"
26 #include "STEER/AliMagFMaps.h"
27 #include "STRUCT/AliBODY.h"
28 #include "STRUCT/AliMAG.h"
29 #include "STRUCT/AliABSOv0.h"
30 #include "STRUCT/AliDIPOv2.h"
31 #include "STRUCT/AliHALL.h"
32 #include "STRUCT/AliFRAMEv2.h"
33 #include "STRUCT/AliSHILv2.h"
34 #include "STRUCT/AliPIPEv0.h"
35 #include "ITS/AliITSvPPRasymmFMD.h"
36 #include "TPC/AliTPCv2.h"
37 #include "TOF/AliTOFv4T0.h"
38 #include "RICH/AliRICHv1.h"
39 #include "ZDC/AliZDCv2.h"
40 #include "TRD/AliTRDv1.h"
41 #include "FMD/AliFMDv1.h"
42 #include "MUON/AliMUONv1.h"
43 #include "MUON/AliMUONSt1GeometryBuilder.h"
44 #include "MUON/AliMUONSt2GeometryBuilder.h"
45 #include "MUON/AliMUONSlatGeometryBuilder.h"
46 #include "MUON/AliMUONTriggerGeometryBuilder.h"
47 #include "PHOS/AliPHOSv1.h"
48 #include "PMD/AliPMDv1.h"
49 #include "START/AliSTARTv1.h"
50 #include "EMCAL/AliEMCALv1.h"
51 #include "CRT/AliCRTv0.h"
52 #include "VZERO/AliVZEROv2.h"
53 #endif
54
55 enum PprRun_t 
56 {
57     test50,
58     kParam_8000,   kParam_4000,  kParam_2000, 
59     kHijing_cent1, kHijing_cent2, 
60     kHijing_per1,  kHijing_per2, kHijing_per3, kHijing_per4,  kHijing_per5,
61     kHijing_jj25,  kHijing_jj50, kHijing_jj75, kHijing_jj100, kHijing_jj200, 
62     kHijing_gj25,  kHijing_gj50, kHijing_gj75, kHijing_gj100, kHijing_gj200,
63     kHijing_pA, kPythia6, 
64     kPythia6Jets20_24,   kPythia6Jets24_29,   kPythia6Jets29_35,
65     kPythia6Jets35_42,   kPythia6Jets42_50,   kPythia6Jets50_60,
66     kPythia6Jets60_72,   kPythia6Jets72_86,   kPythia6Jets86_104,
67     kPythia6Jets104_125, kPythia6Jets125_150, kPythia6Jets150_180,
68     kD0PbPb5500, kD_TRD, kB_TRD, kJpsi_TRD,
69     kU_TRD, kPyJJ, kPyGJ, kRunMax
70 };
71
72 const char* pprRunName[kRunMax] = {
73     "test50",
74     "kParam_8000",   "kParam_4000",  "kParam_2000", 
75     "kHijing_cent1", "kHijing_cent2", 
76     "kHijing_per1",  "kHijing_per2", "kHijing_per3", "kHijing_per4",  
77     "kHijing_per5",
78     "kHijing_jj25",  "kHijing_jj50", "kHijing_jj75", "kHijing_jj100", 
79     "kHijing_jj200", 
80     "kHijing_gj25",  "kHijing_gj50", "kHijing_gj75", "kHijing_gj100", 
81     "kHijing_gj200", "kHijing_pA", "kPythia6", 
82     "kPythia6Jets20_24",   "kPythia6Jets24_29",   "kPythia6Jets29_35",
83     "kPythia6Jets35_42",   "kPythia6Jets42_50",   "kPythia6Jets50_60",
84     "kPythia6Jets60_72",   "kPythia6Jets72_86",   "kPythia6Jets86_104",
85     "kPythia6Jets104_125", "kPythia6Jets125_150", "kPythia6Jets150_180",
86      "kD0PbPb5500", "kD_TRD", 
87     "kB_TRD", "kJpsi_TRD",
88     "kU_TRD", "kPyJJ", "kPyGJ"
89 };
90
91 enum PprGeo_t 
92 {
93     kHoles, kNoHoles
94 };
95
96 enum PprRad_t
97 {
98     kGluonRadiation, kNoGluonRadiation
99 };
100
101 enum PprMag_t
102 {
103     k2kG, k4kG, k5kG
104 };
105
106
107 // This part for configuration    
108 //static PprRun_t srun = test50;
109 static PprRun_t srun = kPythia6;
110 static PprGeo_t sgeo = kHoles;
111 static PprRad_t srad = kGluonRadiation;
112 static PprMag_t smag = k5kG;
113 static Int_t    sseed = 12345; //Set 0 to use the current time
114
115 // Comment line 
116 static TString  comment;
117
118 // Functions
119 Float_t EtaToTheta(Float_t arg);
120 AliGenerator* GeneratorFactory(PprRun_t srun);
121 AliGenHijing* HijingStandard();
122 void ProcessEnvironmentVars();
123
124 void Config()
125 {
126     // ThetaRange is (0., 180.). It was (0.28,179.72) 7/12/00 09:00
127     // Theta range given through pseudorapidity limits 22/6/2001
128
129     // Get settings from environment variables
130     ProcessEnvironmentVars();
131
132     // Set Random Number seed
133     gRandom->SetSeed(sseed);
134     cout<<"Seed for random number generation= "<<gRandom->GetSeed()<<endl; 
135
136
137    // libraries required by geant321
138 #if defined(__CINT__)
139     gSystem->Load("libgeant321");
140 #endif
141
142     new     TGeant3("C++ Interface to Geant3");
143
144     AliRunLoader* rl=0x0;
145
146     cout<<"Config.C: Creating Run Loader ..."<<endl;
147     rl = AliRunLoader::Open("galice.root",
148                             AliConfig::fgkDefaultEventFolderName,
149                             "recreate");
150     if (rl == 0x0)
151       {
152         gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
153         return;
154       }
155     rl->SetCompressionLevel(2);
156     rl->SetNumberOfEventsPerFile(3);
157     gAlice->SetRunLoader(rl);
158
159     //
160     // Set External decayer
161     AliDecayer *decayer = new AliDecayerPythia();
162
163     decayer->SetForceDecay(kAll);
164     decayer->Init();
165     gMC->SetExternalDecayer(decayer);
166     //
167     //
168     //=======================================================================
169     //
170     //=======================================================================
171     // ************* STEERING parameters FOR ALICE SIMULATION **************
172     // --- Specify event type to be tracked through the ALICE setup
173     // --- All positions are in cm, angles in degrees, and P and E in GeV
174
175     gMC->SetProcess("DCAY",1);
176     gMC->SetProcess("PAIR",1);
177     gMC->SetProcess("COMP",1);
178     gMC->SetProcess("PHOT",1);
179     gMC->SetProcess("PFIS",0);
180     gMC->SetProcess("DRAY",0);
181     gMC->SetProcess("ANNI",1);
182     gMC->SetProcess("BREM",1);
183     gMC->SetProcess("MUNU",1);
184     gMC->SetProcess("CKOV",1);
185     gMC->SetProcess("HADR",1);
186     gMC->SetProcess("LOSS",2);
187     gMC->SetProcess("MULS",1);
188     gMC->SetProcess("RAYL",1);
189
190     Float_t cut = 1.e-3;        // 1MeV cut by default
191     Float_t tofmax = 1.e10;
192
193     gMC->SetCut("CUTGAM", cut);
194     gMC->SetCut("CUTELE", cut);
195     gMC->SetCut("CUTNEU", cut);
196     gMC->SetCut("CUTHAD", cut);
197     gMC->SetCut("CUTMUO", cut);
198     gMC->SetCut("BCUTE",  cut); 
199     gMC->SetCut("BCUTM",  cut); 
200     gMC->SetCut("DCUTE",  cut); 
201     gMC->SetCut("DCUTM",  cut); 
202     gMC->SetCut("PPCUTM", cut);
203     gMC->SetCut("TOFMAX", tofmax); 
204
205
206 // Generator Configuration
207     gAlice->SetDebug(0);
208     AliGenerator* gener = GeneratorFactory(srun);
209     gener->SetOrigin(0, 0, 0);    // vertex position
210     gener->SetSigma(0, 0, 5.3);   // Sigma in (X,Y,Z) (cm) on IP position
211     gener->SetCutVertexZ(1.);     // Truncate at 1 sigma
212     gener->SetVertexSmear(kPerEvent); 
213     gener->SetTrackingFlag(1);
214     gener->Init();
215     
216     if (smag == k2kG) {
217         comment = comment.Append(" | L3 field 0.2 T");
218     } else if (smag == k4kG) {
219         comment = comment.Append(" | L3 field 0.4 T");
220     } else if (smag == k5kG) {
221         comment = comment.Append(" | L3 field 0.5 T");
222     }
223     
224     
225     if (srad == kGluonRadiation)
226     {
227         comment = comment.Append(" | Gluon Radiation On");
228         
229     } else {
230         comment = comment.Append(" | Gluon Radiation Off");
231     }
232
233     if (sgeo == kHoles)
234     {
235         comment = comment.Append(" | Holes for PHOS/RICH");
236         
237     } else {
238         comment = comment.Append(" | No holes for PHOS/RICH");
239     }
240
241     printf("\n \n Comment: %s \n \n", comment.Data());
242     
243     
244 // Field (L3 0.4 T)
245     AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., smag);
246     field->SetL3ConstField(0); //Using const. field in the barrel
247     rl->CdGAFile();
248     gAlice->SetField(field);    
249 //
250     Int_t   iABSO   = 1;
251     Int_t   iDIPO   = 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   iRICH   = 1;
262     Int_t   iSHIL   = 1;
263     Int_t   iSTART  = 1;
264     Int_t   iTOF    = 1;
265     Int_t   iTPC    = 1;
266     Int_t   iTRD    = 1;
267     Int_t   iZDC    = 1;
268     Int_t   iEMCAL  = 1;
269     Int_t   iVZERO  = 1;
270     Int_t   iCRT    = 0;
271
272     //=================== Alice BODY parameters =============================
273     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
274
275
276     if (iMAG)
277     {
278         //=================== MAG parameters ============================
279         // --- Start with Magnet since detector layouts may be depending ---
280         // --- on the selected Magnet dimensions ---
281         AliMAG *MAG = new AliMAG("MAG", "Magnet");
282     }
283
284
285     if (iABSO)
286     {
287         //=================== ABSO parameters ============================
288         AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
289     }
290
291     if (iDIPO)
292     {
293         //=================== DIPO parameters ============================
294
295         AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
296     }
297
298     if (iHALL)
299     {
300         //=================== HALL parameters ============================
301
302         AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
303     }
304
305
306     if (iFRAME)
307     {
308         //=================== FRAME parameters ============================
309
310         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
311         if (sgeo == kHoles) {
312             FRAME->SetHoles(1);
313         } else {
314             FRAME->SetHoles(0);
315         }
316     }
317
318     if (iSHIL)
319     {
320         //=================== SHIL parameters ============================
321
322         AliSHIL *SHIL = new AliSHILv2("SHIL", "Shielding Version 2");
323     }
324
325
326     if (iPIPE)
327     {
328         //=================== PIPE parameters ============================
329
330         AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
331     }
332  
333     if(iITS) {
334
335     //=================== ITS parameters ============================
336     //
337     // As the innermost detector in ALICE, the Inner Tracking System "impacts" on
338     // almost all other detectors. This involves the fact that the ITS geometry
339     // still has several options to be followed in parallel in order to determine
340     // the best set-up which minimizes the induced background. All the geometries
341     // available to date are described in the following. Read carefully the comments
342     // and use the default version (the only one uncommented) unless you are making
343     // comparisons and you know what you are doing. In this case just uncomment the
344     // ITS geometry you want to use and run Aliroot.
345     //
346     // Detailed geometries:         
347     //
348     //
349     //AliITS *ITS  = new AliITSv5symm("ITS","Updated ITS TDR detailed version with symmetric services");
350     //
351     //AliITS *ITS  = new AliITSv5asymm("ITS","Updates ITS TDR detailed version with asymmetric services");
352     //
353         AliITSvPPRasymmFMD *ITS  = new AliITSvPPRasymmFMD("ITS","New ITS PPR detailed version with asymmetric services");
354         ITS->SetMinorVersion(2);  // don't touch this parameter if you're not an ITS developer
355         ITS->SetReadDet(kTRUE);   // don't touch this parameter if you're not an ITS developer
356     //    ITS->SetWriteDet("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymm2.det");  // don't touch this parameter if you're not an ITS developer
357         ITS->SetThicknessDet1(200.);   // detector thickness on layer 1 must be in the range [100,300]
358         ITS->SetThicknessDet2(200.);   // detector thickness on layer 2 must be in the range [100,300]
359         ITS->SetThicknessChip1(200.);  // chip thickness on layer 1 must be in the range [150,300]
360         ITS->SetThicknessChip2(200.);  // chip thickness on layer 2 must be in the range [150,300]
361         ITS->SetRails(0);            // 1 --> rails in ; 0 --> rails out
362         ITS->SetCoolingFluid(1);   // 1 --> water ; 0 --> freon
363
364     // Coarse geometries (warning: no hits are produced with these coarse geometries and they unuseful 
365     // for reconstruction !):
366     //                                                     
367     //
368     //AliITSvPPRcoarseasymm *ITS  = new AliITSvPPRcoarseasymm("ITS","New ITS PPR coarse version with asymmetric services");
369     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
370     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
371     //
372     //AliITS *ITS  = new AliITSvPPRcoarsesymm("ITS","New ITS PPR coarse version with symmetric services");
373     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
374     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
375     //                      
376     //
377     //
378     // Geant3 <-> EUCLID conversion
379     // ============================
380     //
381     // SetEUCLID is a flag to output (=1) or not to output (=0) both geometry and
382     // media to two ASCII files (called by default ITSgeometry.euc and
383     // ITSgeometry.tme) in a format understandable to the CAD system EUCLID.
384     // The default (=0) means that you dont want to use this facility.
385     //
386         ITS->SetEUCLID(0);  
387     }
388
389     if (iTPC)
390     {
391         //============================ TPC parameters ================================
392         // --- This allows the user to specify sectors for the SLOW (TPC geometry 2)
393         // --- Simulator. SecAL (SecAU) <0 means that ALL lower (upper)
394         // --- sectors are specified, any value other than that requires at least one 
395         // --- sector (lower or upper)to be specified!
396         // --- Reminder: sectors 1-24 are lower sectors (1-12 -> z>0, 13-24 -> z<0)
397         // ---           sectors 25-72 are the upper ones (25-48 -> z>0, 49-72 -> z<0)
398         // --- SecLows - number of lower sectors specified (up to 6)
399         // --- SecUps - number of upper sectors specified (up to 12)
400         // --- Sens - sensitive strips for the Slow Simulator !!!
401         // --- This does NOT work if all S or L-sectors are specified, i.e.
402         // --- if SecAL or SecAU < 0
403         //
404         //
405         //-----------------------------------------------------------------------------
406
407         //  gROOT->LoadMacro("SetTPCParam.C");
408         //  AliTPCParam *param = SetTPCParam();
409         AliTPC *TPC = new AliTPCv2("TPC", "Default");
410
411         // All sectors included 
412         TPC->SetSecAL(-1);
413         TPC->SetSecAU(-1);
414
415     }
416
417
418     if (iTOF) {
419         //=================== TOF parameters ============================
420         AliTOF *TOF = new AliTOFv4T0("TOF", "normal TOF");
421     }
422
423
424     if (iRICH)
425     {
426         //=================== RICH parameters ===========================
427         AliRICH *RICH = new AliRICHv1("RICH", "normal RICH");
428
429     }
430
431
432     if (iZDC)
433     {
434         //=================== ZDC parameters ============================
435
436         AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
437     }
438
439     if (iTRD)
440     {
441         //=================== TRD parameters ============================
442
443         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
444
445         // Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 90% Xe + 10% CO2)
446         TRD->SetGasMix(1);
447         if (sgeo == kHoles) {
448             // With hole in front of PHOS
449             TRD->SetPHOShole();
450             // With hole in front of RICH
451             TRD->SetRICHhole();
452         }
453             // Switch on TR
454             AliTRDsim *TRDsim = TRD->CreateTR();
455     }
456
457     if (iFMD)
458     {
459         //=================== FMD parameters ============================
460         AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
461    }
462
463     if (iMUON)
464     {
465         //=================== MUON parameters ===========================
466
467         AliMUON *MUON = new AliMUONv1("MUON", "default");
468         MUON->AddGeometryBuilder(new AliMUONSt1GeometryBuilder(MUON));
469         MUON->AddGeometryBuilder(new AliMUONSt2GeometryBuilder(MUON));
470         MUON->AddGeometryBuilder(new AliMUONSlatGeometryBuilder(MUON));
471         MUON->AddGeometryBuilder(new AliMUONTriggerGeometryBuilder(MUON));
472     }
473     //=================== PHOS parameters ===========================
474
475     if (iPHOS)
476     {
477         AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
478     }
479
480
481     if (iPMD)
482     {
483         //=================== PMD parameters ============================
484         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
485     }
486
487     if (iSTART)
488     {
489         //=================== START parameters ============================
490         AliSTART *START = new AliSTARTv1("START", "START Detector");
491     }
492
493     if (iEMCAL)
494     {
495         //=================== EMCAL parameters ============================
496         AliEMCAL *EMCAL = new AliEMCALv1("EMCAL", "EMCAL_55_25");
497     }
498
499      if (iCRT)
500     {
501         //=================== CRT parameters ============================
502         AliCRT *CRT = new AliCRTv0("CRT", "normal ACORDE");
503     }
504
505      if (iVZERO)
506     {
507         //=================== CRT parameters ============================
508         AliVZERO *VZERO = new AliVZEROv2("VZERO", "normal VZERO");
509     }
510  
511              
512 }
513
514 Float_t EtaToTheta(Float_t arg){
515   return (180./TMath::Pi())*2.*atan(exp(-arg));
516 }
517
518
519
520 AliGenerator* GeneratorFactory(PprRun_t srun) {
521     Int_t isw = 3;
522     if (srad == kNoGluonRadiation) isw = 0;
523     
524
525     AliGenerator * gGener = 0x0;
526     switch (srun) {
527     case test50:
528       {
529         comment = comment.Append(":HIJINGparam test 50 particles");
530         AliGenHIJINGpara *gener = new AliGenHIJINGpara(50);
531         gener->SetMomentumRange(0, 999999.);
532         gener->SetPhiRange(0., 360.);
533         // Set pseudorapidity range from -8 to 8.
534         Float_t thmin = EtaToTheta(8);   // theta min. <---> eta max
535         Float_t thmax = EtaToTheta(-8);  // theta max. <---> eta min 
536         gener->SetThetaRange(thmin,thmax);
537         gGener=gener;
538       }
539       break;
540     case kParam_8000:
541       {
542         comment = comment.Append(":HIJINGparam N=8000");
543         AliGenHIJINGpara *gener = new AliGenHIJINGpara(86030);
544         gener->SetMomentumRange(0, 999999.);
545         gener->SetPhiRange(0., 360.);
546         // Set pseudorapidity range from -8 to 8.
547         Float_t thmin = EtaToTheta(8);   // theta min. <---> eta max
548         Float_t thmax = EtaToTheta(-8);  // theta max. <---> eta min 
549         gener->SetThetaRange(thmin,thmax);
550         gGener=gener;
551       }
552       break;
553     case kParam_4000:
554       {
555         comment = comment.Append("HIJINGparam N=4000");
556         AliGenHIJINGpara *gener = new AliGenHIJINGpara(43015);
557         gener->SetMomentumRange(0, 999999.);
558         gener->SetPhiRange(0., 360.);
559         // Set pseudorapidity range from -8 to 8.
560         Float_t thmin = EtaToTheta(8);   // theta min. <---> eta max
561         Float_t thmax = EtaToTheta(-8);  // theta max. <---> eta min 
562         gener->SetThetaRange(thmin,thmax);
563         gGener=gener;
564       }
565         break;
566     case kParam_2000:
567       {
568         comment = comment.Append("HIJINGparam N=2000");
569         AliGenHIJINGpara *gener = new AliGenHIJINGpara(21507);
570         gener->SetMomentumRange(0, 999999.);
571         gener->SetPhiRange(0., 360.);
572         // Set pseudorapidity range from -8 to 8.
573         Float_t thmin = EtaToTheta(8);   // theta min. <---> eta max
574         Float_t thmax = EtaToTheta(-8);  // theta max. <---> eta min 
575         gener->SetThetaRange(thmin,thmax);
576         gGener=gener;
577       }
578       break;
579 //
580 //  Hijing Central
581 //
582     case kHijing_cent1:
583       {
584         comment = comment.Append("HIJING cent1");
585         AliGenHijing *gener = HijingStandard();
586 // impact parameter range
587         gener->SetImpactParameterRange(0., 5.);
588         gGener=gener;
589       }
590       break;
591     case kHijing_cent2:
592       {
593         comment = comment.Append("HIJING cent2");
594         AliGenHijing *gener = HijingStandard();
595 // impact parameter range
596         gener->SetImpactParameterRange(0., 2.);
597         gGener=gener;
598       }
599       break;
600 //
601 // Hijing Peripheral 
602 //
603     case kHijing_per1:
604       {
605         comment = comment.Append("HIJING per1");
606         AliGenHijing *gener = HijingStandard();
607 // impact parameter range
608         gener->SetImpactParameterRange(5., 8.6);
609         gGener=gener;
610       }
611       break;
612     case kHijing_per2:
613       {
614         comment = comment.Append("HIJING per2");
615         AliGenHijing *gener = HijingStandard();
616 // impact parameter range
617         gener->SetImpactParameterRange(8.6, 11.2);
618         gGener=gener;
619       }
620       break;
621     case kHijing_per3:
622       {
623         comment = comment.Append("HIJING per3");
624         AliGenHijing *gener = HijingStandard();
625 // impact parameter range
626         gener->SetImpactParameterRange(11.2, 13.2);
627         gGener=gener;
628       }
629       break;
630     case kHijing_per4:
631       {
632         comment = comment.Append("HIJING per4");
633         AliGenHijing *gener = HijingStandard();
634 // impact parameter range
635         gener->SetImpactParameterRange(13.2, 15.);
636         gGener=gener;
637       }
638       break;
639     case kHijing_per5:
640       {
641         comment = comment.Append("HIJING per5");
642         AliGenHijing *gener = HijingStandard();
643 // impact parameter range
644         gener->SetImpactParameterRange(15., 100.);
645         gGener=gener;
646       }
647       break;
648 //
649 //  Jet-Jet
650 //
651     case kHijing_jj25:
652       {
653         comment = comment.Append("HIJING Jet 25 GeV");
654         AliGenHijing *gener = HijingStandard();
655 // impact parameter range
656         gener->SetImpactParameterRange(0., 5.);
657         // trigger
658         gener->SetTrigger(1);
659         gener->SetPtJet(25.);
660         gener->SetRadiation(isw);
661         gener->SetSimpleJets(!isw);
662         gener->SetJetEtaRange(-0.3,0.3);
663         gener->SetJetPhiRange(75., 165.);   
664         gGener=gener;
665       }
666       break;
667
668     case kHijing_jj50:
669       {
670         comment = comment.Append("HIJING Jet 50 GeV");
671         AliGenHijing *gener = HijingStandard();
672 // impact parameter range
673         gener->SetImpactParameterRange(0., 5.);
674         // trigger
675         gener->SetTrigger(1);
676         gener->SetPtJet(50.);
677         gener->SetRadiation(isw);
678         gener->SetSimpleJets(!isw);
679         gener->SetJetEtaRange(-0.3,0.3);
680         gener->SetJetPhiRange(75., 165.);   
681         gGener=gener;
682       }
683         break;
684
685     case kHijing_jj75:
686       {
687         comment = comment.Append("HIJING Jet 75 GeV");
688         AliGenHijing *gener = HijingStandard();
689 // impact parameter range
690         gener->SetImpactParameterRange(0., 5.);
691         // trigger
692         gener->SetTrigger(1);
693         gener->SetPtJet(75.);
694         gener->SetRadiation(isw);
695         gener->SetSimpleJets(!isw);
696         gener->SetJetEtaRange(-0.3,0.3);
697         gener->SetJetPhiRange(75., 165.);   
698         gGener=gener;
699       }
700       break;
701
702     case kHijing_jj100:
703       {
704         comment = comment.Append("HIJING Jet 100 GeV");
705         AliGenHijing *gener = HijingStandard();
706 // impact parameter range
707         gener->SetImpactParameterRange(0., 5.);
708         // trigger
709         gener->SetTrigger(1);
710         gener->SetPtJet(100.);
711         gener->SetRadiation(isw);
712         gener->SetSimpleJets(!isw);
713         gener->SetJetEtaRange(-0.3,0.3);
714         gener->SetJetPhiRange(75., 165.);   
715         gGener=gener;
716       }
717       break;
718
719     case kHijing_jj200:
720       {
721         comment = comment.Append("HIJING Jet 200 GeV");
722         AliGenHijing *gener = HijingStandard();
723 // impact parameter range
724         gener->SetImpactParameterRange(0., 5.);
725         // trigger
726         gener->SetTrigger(1);
727         gener->SetPtJet(200.);
728         gener->SetRadiation(isw);
729         gener->SetSimpleJets(!isw);
730         gener->SetJetEtaRange(-0.3,0.3);
731         gener->SetJetPhiRange(75., 165.);   
732         gGener=gener;
733       }
734       break;
735 //
736 // Gamma-Jet
737 //
738     case kHijing_gj25:
739       {
740         comment = comment.Append("HIJING Gamma 25 GeV");
741         AliGenHijing *gener = HijingStandard();
742 // impact parameter range
743         gener->SetImpactParameterRange(0., 5.);
744         // trigger
745         gener->SetTrigger(2);
746         gener->SetPtJet(25.);
747         gener->SetRadiation(isw);
748         gener->SetSimpleJets(!isw);
749         gener->SetJetEtaRange(-0.12, 0.12);
750         gener->SetJetPhiRange(220., 320.);
751         gGener=gener;
752       }
753       break;
754
755     case kHijing_gj50:
756       {
757         comment = comment.Append("HIJING Gamma 50 GeV");
758         AliGenHijing *gener = HijingStandard();
759 // impact parameter range
760         gener->SetImpactParameterRange(0., 5.);
761         // trigger
762         gener->SetTrigger(2);
763         gener->SetPtJet(50.);
764         gener->SetRadiation(isw);
765         gener->SetSimpleJets(!isw);
766         gener->SetJetEtaRange(-0.12, 0.12);
767         gener->SetJetPhiRange(220., 320.);
768         gGener=gener;
769       }
770       break;
771
772     case kHijing_gj75:
773       {
774         comment = comment.Append("HIJING Gamma 75 GeV");
775         AliGenHijing *gener = HijingStandard();
776 // impact parameter range
777         gener->SetImpactParameterRange(0., 5.);
778         // trigger
779         gener->SetTrigger(2);
780         gener->SetPtJet(75.);
781         gener->SetRadiation(isw);
782         gener->SetSimpleJets(!isw);
783         gener->SetJetEtaRange(-0.12, 0.12);
784         gener->SetJetPhiRange(220., 320.);
785         gGener=gener;
786       }
787       break;
788
789     case kHijing_gj100:
790       {
791         comment = comment.Append("HIJING Gamma 100 GeV");
792         AliGenHijing *gener = HijingStandard();
793 // impact parameter range
794         gener->SetImpactParameterRange(0., 5.);
795         // trigger
796         gener->SetTrigger(2);
797         gener->SetPtJet(100.);
798         gener->SetRadiation(isw);
799         gener->SetSimpleJets(!isw);
800         gener->SetJetEtaRange(-0.12, 0.12);
801         gener->SetJetPhiRange(220., 320.);
802         gGener=gener;
803       }
804       break;
805
806     case kHijing_gj200:
807       {
808         comment = comment.Append("HIJING Gamma 200 GeV");
809         AliGenHijing *gener = HijingStandard();
810 // impact parameter range
811         gener->SetImpactParameterRange(0., 5.);
812         // trigger
813         gener->SetTrigger(2);
814         gener->SetPtJet(200.);
815         gener->SetRadiation(isw);
816         gener->SetSimpleJets(!isw);
817         gener->SetJetEtaRange(-0.12, 0.12);
818         gener->SetJetPhiRange(220., 320.);
819         gGener=gener;
820       }
821       break;
822     case kHijing_pA:
823       {
824         comment = comment.Append("HIJING pA");
825
826         AliGenCocktail *gener  = new AliGenCocktail();
827
828         AliGenHijing   *hijing = new AliGenHijing(-1);
829 // centre of mass energy 
830         hijing->SetEnergyCMS(TMath::Sqrt(82./208.) * 14000.);
831 // impact parameter range
832         hijing->SetImpactParameterRange(0., 15.);
833 // reference frame
834         hijing->SetReferenceFrame("CMS");
835         hijing->SetBoostLHC(1);
836 // projectile
837         hijing->SetProjectile("P", 1, 1);
838         hijing->SetTarget    ("A", 208, 82);
839 // tell hijing to keep the full parent child chain
840         hijing->KeepFullEvent();
841 // enable jet quenching
842         hijing->SetJetQuenching(0);
843 // enable shadowing
844         hijing->SetShadowing(1);
845 // Don't track spectators
846         hijing->SetSpectators(0);
847 // kinematic selection
848         hijing->SetSelectAll(0);
849 //
850         AliGenSlowNucleons*  gray    = new AliGenSlowNucleons(1);
851         AliSlowNucleonModel* model   = new AliSlowNucleonModelExp();
852         gray->SetSlowNucleonModel(model);
853         gray->SetDebug(1);
854         gener->AddGenerator(hijing,"Hijing pPb", 1);
855         gener->AddGenerator(gray,  "Gray Particles",1);
856         gGener=gener;
857       }
858       break;
859     case kPythia6:
860       {
861         comment = comment.Append(":Pythia p-p @ 14 TeV");
862         AliGenPythia *gener = new AliGenPythia(-1); 
863         gener->SetMomentumRange(0,999999);
864         gener->SetPhiRange(-180,180);
865         gener->SetThetaRange(0., 180.);
866         gener->SetYRange(-12,12);
867         gener->SetPtRange(0,1000);
868         gener->SetStrucFunc(kCTEQ4L);   
869         gener->SetProcess(kPyMb);
870         gener->SetEnergyCMS(14000.);
871         gGener=gener;
872       }
873       break;
874     case kPythia6Jets20_24:
875       {
876         comment = comment.Append(":Pythia jets 20-24 GeV @ 5.5 TeV");
877         AliGenPythia * gener = new AliGenPythia(-1);
878         gener->SetEnergyCMS(5500.);//        Centre of mass energy
879         gener->SetProcess(kPyJets);//        Process type
880         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
881         gener->SetJetPhiRange(0., 360.);
882         gener->SetJetEtRange(10., 1000.);
883         gener->SetGluonRadiation(1,1);
884         //    gener->SetPtKick(0.);
885         //   Structure function
886         gener->SetStrucFunc(kCTEQ4L);
887         gener->SetPtHard(20., 24.);// Pt transfer of the hard scattering
888         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
889         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
890         gGener=gener;
891       }
892       break;
893     case kPythia6Jets24_29:
894       {
895         comment = comment.Append(":Pythia jets 24-29 GeV @ 5.5 TeV");
896         AliGenPythia * gener = new AliGenPythia(-1);
897         gener->SetEnergyCMS(5500.);//        Centre of mass energy
898         gener->SetProcess(kPyJets);//        Process type
899         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
900         gener->SetJetPhiRange(0., 360.);
901         gener->SetJetEtRange(10., 1000.);
902         gener->SetGluonRadiation(1,1);
903         //    gener->SetPtKick(0.);
904         //   Structure function
905         gener->SetStrucFunc(kCTEQ4L);
906         gener->SetPtHard(24., 29.);// Pt transfer of the hard scattering
907         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
908         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
909         gGener=gener;
910       }
911       break;
912     case kPythia6Jets29_35:
913       {
914         comment = comment.Append(":Pythia jets 29-35 GeV @ 5.5 TeV");
915         AliGenPythia * gener = new AliGenPythia(-1);
916         gener->SetEnergyCMS(5500.);//        Centre of mass energy
917         gener->SetProcess(kPyJets);//        Process type
918         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
919         gener->SetJetPhiRange(0., 360.);
920         gener->SetJetEtRange(10., 1000.);
921         gener->SetGluonRadiation(1,1);
922         //    gener->SetPtKick(0.);
923         //   Structure function
924         gener->SetStrucFunc(kCTEQ4L);
925         gener->SetPtHard(29., 35.);// Pt transfer of the hard scattering
926         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
927         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
928         gGener=gener;
929       }
930       break;
931     case kPythia6Jets35_42:
932       {
933         comment = comment.Append(":Pythia jets 35-42 GeV @ 5.5 TeV");
934         AliGenPythia * gener = new AliGenPythia(-1);
935         gener->SetEnergyCMS(5500.);//        Centre of mass energy
936         gener->SetProcess(kPyJets);//        Process type
937         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
938         gener->SetJetPhiRange(0., 360.);
939         gener->SetJetEtRange(10., 1000.);
940         gener->SetGluonRadiation(1,1);
941         //    gener->SetPtKick(0.);
942         //   Structure function
943         gener->SetStrucFunc(kCTEQ4L);
944         gener->SetPtHard(35., 42.);// Pt transfer of the hard scattering
945         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
946         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
947         gGener=gener;
948       }
949       break;
950     case kPythia6Jets42_50:
951       {
952         comment = comment.Append(":Pythia jets 42-50 GeV @ 5.5 TeV");
953         AliGenPythia * gener = new AliGenPythia(-1);
954         gener->SetEnergyCMS(5500.);//        Centre of mass energy
955         gener->SetProcess(kPyJets);//        Process type
956         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
957         gener->SetJetPhiRange(0., 360.);
958         gener->SetJetEtRange(10., 1000.);
959         gener->SetGluonRadiation(1,1);
960         //    gener->SetPtKick(0.);
961         //   Structure function
962         gener->SetStrucFunc(kCTEQ4L);
963         gener->SetPtHard(42., 50.);// Pt transfer of the hard scattering
964         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
965         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
966         gGener=gener;
967       }
968       break;
969     case kPythia6Jets50_60:
970       {
971         comment = comment.Append(":Pythia jets 50-60 GeV @ 5.5 TeV");
972         AliGenPythia * gener = new AliGenPythia(-1);
973         gener->SetEnergyCMS(5500.);//        Centre of mass energy
974         gener->SetProcess(kPyJets);//        Process type
975         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
976         gener->SetJetPhiRange(0., 360.);
977         gener->SetJetEtRange(10., 1000.);
978         gener->SetGluonRadiation(1,1);
979         //    gener->SetPtKick(0.);
980         //   Structure function
981         gener->SetStrucFunc(kCTEQ4L);
982         gener->SetPtHard(50., 60.);// Pt transfer of the hard scattering
983         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
984         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
985         gGener=gener;
986       }
987       break;
988     case kPythia6Jets60_72:
989       {
990         comment = comment.Append(":Pythia jets 60-72 GeV @ 5.5 TeV");
991         AliGenPythia * gener = new AliGenPythia(-1);
992         gener->SetEnergyCMS(5500.);//        Centre of mass energy
993         gener->SetProcess(kPyJets);//        Process type
994         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
995         gener->SetJetPhiRange(0., 360.);
996         gener->SetJetEtRange(10., 1000.);
997         gener->SetGluonRadiation(1,1);
998         //    gener->SetPtKick(0.);
999         //   Structure function
1000         gener->SetStrucFunc(kCTEQ4L);
1001         gener->SetPtHard(60., 72.);// Pt transfer of the hard scattering
1002         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1003         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1004         gGener=gener;
1005       }
1006       break;
1007     case kPythia6Jets72_86:
1008       {
1009         comment = comment.Append(":Pythia jets 72-86 GeV @ 5.5 TeV");
1010         AliGenPythia * gener = new AliGenPythia(-1);
1011         gener->SetEnergyCMS(5500.);//        Centre of mass energy
1012         gener->SetProcess(kPyJets);//        Process type
1013         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
1014         gener->SetJetPhiRange(0., 360.);
1015         gener->SetJetEtRange(10., 1000.);
1016         gener->SetGluonRadiation(1,1);
1017         //    gener->SetPtKick(0.);
1018         //   Structure function
1019         gener->SetStrucFunc(kCTEQ4L);
1020         gener->SetPtHard(72., 86.);// Pt transfer of the hard scattering
1021         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1022         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1023         gGener=gener;
1024       }
1025       break;
1026     case kPythia6Jets86_104:
1027       {
1028         comment = comment.Append(":Pythia jets 86-104 GeV @ 5.5 TeV");
1029         AliGenPythia * gener = new AliGenPythia(-1);
1030         gener->SetEnergyCMS(5500.);//        Centre of mass energy
1031         gener->SetProcess(kPyJets);//        Process type
1032         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
1033         gener->SetJetPhiRange(0., 360.);
1034         gener->SetJetEtRange(10., 1000.);
1035         gener->SetGluonRadiation(1,1);
1036         //    gener->SetPtKick(0.);
1037         //   Structure function
1038         gener->SetStrucFunc(kCTEQ4L);
1039         gener->SetPtHard(86., 104.);// Pt transfer of the hard scattering
1040         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1041         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1042         gGener=gener;
1043       }
1044       break;
1045     case kPythia6Jets104_125:
1046       {
1047         comment = comment.Append(":Pythia jets 105-125 GeV @ 5.5 TeV");
1048         AliGenPythia * gener = new AliGenPythia(-1);
1049         gener->SetEnergyCMS(5500.);//        Centre of mass energy
1050         gener->SetProcess(kPyJets);//        Process type
1051         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
1052         gener->SetJetPhiRange(0., 360.);
1053         gener->SetJetEtRange(10., 1000.);
1054         gener->SetGluonRadiation(1,1);
1055         //    gener->SetPtKick(0.);
1056         //   Structure function
1057         gener->SetStrucFunc(kCTEQ4L);
1058         gener->SetPtHard(104., 125.);// Pt transfer of the hard scattering
1059         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1060         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1061         gGener=gener;
1062       }
1063       break;
1064     case kPythia6Jets125_150:
1065       {
1066         comment = comment.Append(":Pythia jets 125-150 GeV @ 5.5 TeV");
1067         AliGenPythia * gener = new AliGenPythia(-1);
1068         gener->SetEnergyCMS(5500.);//        Centre of mass energy
1069         gener->SetProcess(kPyJets);//        Process type
1070         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
1071         gener->SetJetPhiRange(0., 360.);
1072         gener->SetJetEtRange(10., 1000.);
1073         gener->SetGluonRadiation(1,1);
1074         //    gener->SetPtKick(0.);
1075         //   Structure function
1076         gener->SetStrucFunc(kCTEQ4L);
1077         gener->SetPtHard(125., 150.);// Pt transfer of the hard scattering
1078         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1079         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1080         gGener=gener;
1081       }
1082       break;
1083     case kPythia6Jets150_180:
1084       {
1085         comment = comment.Append(":Pythia jets 150-180 GeV @ 5.5 TeV");
1086         AliGenPythia * gener = new AliGenPythia(-1);
1087         gener->SetEnergyCMS(5500.);//        Centre of mass energy
1088         gener->SetProcess(kPyJets);//        Process type
1089         gener->SetJetEtaRange(-0.5, 0.5);//  Final state kinematic cuts
1090         gener->SetJetPhiRange(0., 360.);
1091         gener->SetJetEtRange(10., 1000.);
1092         gener->SetGluonRadiation(1,1);
1093         //    gener->SetPtKick(0.);
1094         //   Structure function
1095         gener->SetStrucFunc(kCTEQ4L);
1096         gener->SetPtHard(150., 180.);// Pt transfer of the hard scattering
1097         gener->SetPycellParameters(2., 274, 432, 0., 4., 5., 1.0);
1098         gener->SetForceDecay(kAll);//  Decay type (semielectronic, etc.)
1099         gGener=gener;
1100       }
1101       break;
1102     case kD0PbPb5500:
1103       {
1104         comment = comment.Append(" D0 in Pb-Pb at 5.5 TeV");
1105         AliGenPythia * gener = new AliGenPythia(10);
1106         gener->SetProcess(kPyD0PbPbMNR);
1107         gener->SetStrucFunc(kCTEQ4L);
1108         gener->SetPtHard(2.1,-1.0);
1109         gener->SetEnergyCMS(5500.);
1110         gener->SetNuclei(208,208);
1111         gGener=gener;
1112       }
1113       break;
1114     case kD_TRD:
1115       {
1116         comment = comment.Append(" D0 for TRD at 5.5 TeV");
1117         AliGenPythia *gener = new AliGenPythia(1);
1118         gener->SetCutOnChild(0);
1119         gener->SetStrucFunc(kCTEQ4L);
1120         gener->SetProcess(kPyCharm);
1121         gener->SetPtHard(0.,-1);
1122         gener->SetEnergyCMS(5500.);
1123         gGener=gener;
1124       }
1125       break;
1126     case kB_TRD:
1127       {
1128         comment = comment.Append(" B for TRD at 5.5 TeV");
1129         AliGenPythia *gener = new AliGenPythia(1);
1130         gener->SetCutOnChild(0);
1131         gener->SetStrucFunc(kCTEQ4L);
1132         gener->SetProcess(kPyBeauty);
1133         gener->SetPtHard(0.,-1);
1134         gener->SetEnergyCMS(5500.);
1135         gGener=gener;
1136       }
1137       break;
1138     case kJpsi_TRD:
1139       {
1140         comment = comment.Append(" J/psi for TRD at 5.5 TeV");
1141         AliGenParam *gener = new AliGenParam(1,new AliGenGSIlib(),
1142                                              AliGenGSIlib::kJPsi,"MUON");
1143         gener->SetMomentumRange(0,999);
1144         gener->SetPtRange(0,30.);
1145         gener->SetPhiRange(0., 360.);
1146         gener->SetYRange(-0.9,+0.9);
1147         gener->SetForceDecay(kDiElectron);
1148         gGener=gener;
1149       }
1150       break;
1151     case kU_TRD:
1152       {
1153         comment = comment.Append(" Upsilon for TRD at 5.5 TeV");
1154         AliGenParam *gener = new AliGenParam(1,new AliGenGSIlib(),
1155                                              AliGenGSIlib::kUpsilon,"RITMAN");
1156         gener->SetMomentumRange(0,999);
1157         gener->SetPtRange(0,30.);
1158         gener->SetPhiRange(0., 360.);
1159         gener->SetYRange(-0.9,0.9);
1160         gener->SetForceDecay(kDiElectron);
1161         gGener=gener;
1162       }
1163       break;
1164     case kPyJJ:
1165       {
1166         comment = comment.Append(" Jet-jet at 5.5 TeV");
1167         AliGenPythia *gener = new AliGenPythia(-1);
1168         gener->SetEnergyCMS(5500.);
1169         gener->SetProcess(kPyJets);
1170         Double_t ptHardMin=10.0, ptHardMax=-1.0;
1171         gener->SetPtHard(ptHardMin,ptHardMax);
1172         gener->SetYHard(-0.7,0.7);
1173         gener->SetJetEtaRange(-0.2,0.2);
1174         gener->SetEventListRange(0,1);
1175         gGener=gener;
1176       }
1177       break;
1178     case kPyGJ:
1179       {
1180         comment = comment.Append(" Gamma-jet at 5.5 TeV");
1181         AliGenPythia *gener = new AliGenPythia(-1);
1182         gener->SetEnergyCMS(5500.);
1183         gener->SetProcess(kPyDirectGamma);
1184         Double_t ptHardMin=10.0, ptHardMax=-1.0;
1185         gener->SetPtHard(ptHardMin,ptHardMax);
1186         gener->SetYHard(-1.0,1.0);
1187         gener->SetGammaEtaRange(-0.13,0.13);
1188         gener->SetGammaPhiRange(210.,330.);
1189         gener->SetEventListRange(0,1);
1190         gGener=gener;
1191       }
1192       break;
1193     default: break;
1194     }
1195     return gGener;
1196 }
1197
1198 AliGenHijing* HijingStandard()
1199 {
1200     AliGenHijing *gener = new AliGenHijing(-1);
1201 // centre of mass energy 
1202     gener->SetEnergyCMS(5500.);
1203 // reference frame
1204     gener->SetReferenceFrame("CMS");
1205 // projectile
1206      gener->SetProjectile("A", 208, 82);
1207      gener->SetTarget    ("A", 208, 82);
1208 // tell hijing to keep the full parent child chain
1209      gener->KeepFullEvent();
1210 // enable jet quenching
1211      gener->SetJetQuenching(1);
1212 // enable shadowing
1213      gener->SetShadowing(1);
1214 // neutral pion and heavy particle decays switched off
1215      gener->SetDecaysOff(1);
1216 // Don't track spectators
1217      gener->SetSpectators(0);
1218 // kinematic selection
1219      gener->SetSelectAll(0);
1220      return gener;
1221 }
1222
1223
1224 void ProcessEnvironmentVars()
1225 {
1226     // Run type
1227     if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
1228       for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
1229         if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0) {
1230           srun = (PprRun_t)iRun;
1231           cout<<"Run type set to "<<pprRunName[iRun]<<endl;
1232         }
1233       }
1234     }
1235
1236     // Random Number seed
1237     if (gSystem->Getenv("CONFIG_SEED")) {
1238       sseed = atoi(gSystem->Getenv("CONFIG_SEED"));
1239     }
1240 }