update for new geant4 production
[u/mrichter/AliRoot.git] / test / vmctest / production / Config.C
1 // $Id$
2 //
3 // Configuration for the Geant4 production 2010
4 // By E. Sicking, CERN
5 //
6
7 #if !defined(__CINT__) || defined(__MAKECINT__)
8 #include <Riostream.h>
9 #include <TRandom.h>
10 #include <TDatime.h>
11 #include <TSystem.h>
12 #include <TVirtualMC.h>
13 #include <TGeant3TGeo.h>
14 #include "STEER/AliRunLoader.h"
15 #include "STEER/AliRun.h"
16 #include "STEER/AliConfig.h"
17 #include "PYTHIA6/AliDecayerPythia.h"
18 #include "PYTHIA6/AliGenPythia.h"
19 #include "TDPMjet/AliGenDPMjet.h"
20 #include "STEER/AliMagFCheb.h"
21 #include "STRUCT/AliBODY.h"
22 #include "STRUCT/AliMAG.h"
23 #include "STRUCT/AliABSOv3.h"
24 #include "STRUCT/AliDIPOv3.h"
25 #include "STRUCT/AliHALLv3.h"
26 #include "STRUCT/AliFRAMEv2.h"
27 #include "STRUCT/AliSHILv3.h"
28 #include "STRUCT/AliPIPEv3.h"
29 #include "ITS/AliITSv11Hybrid.h"
30 #include "TPC/AliTPCv2.h"
31 #include "TOF/AliTOFv6T0.h"
32 #include "HMPID/AliHMPIDv3.h"
33 #include "ZDC/AliZDCv3.h"
34 #include "TRD/AliTRDv1.h"
35 #include "TRD/AliTRDgeometry.h"
36 #include "FMD/AliFMDv1.h"
37 #include "MUON/AliMUONv1.h"
38 #include "PHOS/AliPHOSv1.h"
39 #include "PHOS/AliPHOSSimParam.h"
40 #include "PMD/AliPMDv1.h"
41 #include "T0/AliT0v1.h"
42 #include "EMCAL/AliEMCALv2.h"
43 #include "ACORDE/AliACORDEv1.h"
44 #include "VZERO/AliVZEROv7.h"
45 #endif
46
47
48 enum PDC06Proc_t 
49 {
50   kPythia6, kPythia6D6T, kPythia6D6T_Flat, kPythia6ATLAS, kPythia6ATLAS_Flat, kPythia6_Perugia0, kPhojet, kRunMax
51 };
52
53 const char * pprRunName[] = {
54     "kPythia6", "kPythia6D6T", "kPythia6D6T_Flat", "kPythia6ATLAS", "kPythia6ATLAS_Flat", "kPythia6_Perugia0", "kPhojet" 
55 };
56
57 enum Mag_t
58 {
59   kNoField, k5kG, kFieldMax
60 };
61
62 const char * pprField[] = {
63   "kNoField", "k5kG"
64 };
65
66 enum PhysicsList_t 
67 {
68   QGSP_BERT_EMV,CHIPS,QGSP_BERT_CHIPS,QGSP_BERT_EMV_OPTICAL, CHIPS_OPTICAL, QGSP_BERT_CHIPS_OPTICAL, kListMax
69 };
70
71 const char * physicsListName[] = {
72   "QGSP_BERT_EMV",         "CHIPS",         "QGSP_BERT_CHIPS", 
73   "QGSP_BERT_EMV_OPTICAL", "CHIPS_OPTICAL", "QGSP_BERT_CHIPS_OPTICAL"
74 };
75
76 enum PprTrigConf_t
77 {
78     kDefaultPPTrig, kDefaultPbPbTrig
79 };
80
81 const char * pprTrigConfName[] = {
82     "p-p","Pb-Pb"
83 };
84
85 //--- Functions ---
86 class AliGenPythia;
87 AliGenerator *MbPythia();
88 AliGenerator *MbPythiaTuneD6T();
89 AliGenerator *MbPhojet();
90 void ProcessEnvironmentVars();
91
92 // Geterator, field, beam energy
93 static PDC06Proc_t   proc         = kPhojet;
94 static Mag_t         mag          = k5kG;
95 static Float_t       energy       = 10000; // energy in CMS
96 static PhysicsList_t physicslist  = QGSP_BERT_EMV;
97 static PprTrigConf_t strig = kDefaultPPTrig; // default pp trigger configuration
98
99 //========================//
100 // Set Random Number seed //
101 //========================//
102 TDatime dt;
103 static UInt_t seed    = dt.Get();
104
105 // Comment line
106 static TString comment;
107
108 void Config()
109 {
110     
111
112   // Get settings from environment variables
113   ProcessEnvironmentVars();
114
115   gRandom->SetSeed(seed);
116   cerr<<"Seed for random number generation= "<<seed<<endl; 
117
118   // Libraries required by geant321
119 #if defined(__CINT__)
120   gSystem->Load("liblhapdf");       // Parton density functions
121   gSystem->Load("libEGPythia6");    // TGenerator interface
122
123   if (proc == kPythia6 || proc == kPhojet) {
124       gSystem->Load("libpythia6");        // Pythia 6.2 
125   } else {
126       gSystem->Load("libpythia6.4.21");   // Pythia 6.4
127   } 
128
129   gSystem->Load("libAliPythia6");   // ALICE specific implementations
130   //gSystem->Load("libgeant321");
131 #endif
132
133   //  new TGeant3TGeo("C++ Interface to Geant3");
134
135   //=======================================================================
136   //  Create the output file
137
138    
139   AliRunLoader* rl=0x0;
140
141   cout<<"Config.C: Creating Run Loader ..."<<endl;
142   rl = AliRunLoader::Open("galice.root",
143                           AliConfig::GetDefaultEventFolderName(),
144                           "recreate");
145   if (rl == 0x0)
146     {
147       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
148       return;
149     }
150   rl->SetCompressionLevel(2);
151   rl->SetNumberOfEventsPerFile(1000);
152   gAlice->SetRunLoader(rl);
153   // gAlice->SetGeometryFromFile("geometry.root");
154   // gAlice->SetGeometryFromCDB();
155   
156   // Set the trigger configuration: proton-proton
157   AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
158   cout<<"Trigger configuration is set to  "<<pprTrigConfName[strig]<<endl;
159
160   printf("\n \n Comment: %s \n \n", comment.Data());
161  
162   rl->CdGAFile();
163   
164   Int_t iABSO  = 1;
165   Int_t iACORDE= 0;
166   Int_t iDIPO  = 1;
167   Int_t iEMCAL = 1;
168   Int_t iFMD   = 1;
169   Int_t iFRAME = 1;
170   Int_t iHALL  = 1;
171   Int_t iITS   = 1;
172   Int_t iMAG   = 1;
173   Int_t iMUON  = 1;
174   Int_t iPHOS  = 1;
175   Int_t iPIPE  = 1;
176   Int_t iPMD   = 1;
177   Int_t iHMPID = 1;
178   Int_t iSHIL  = 1;
179   Int_t iT0    = 1;
180   Int_t iTOF   = 1;
181   Int_t iTPC   = 1;
182   Int_t iTRD   = 1;
183   Int_t iVZERO = 1;
184   Int_t iZDC   = 1;
185   
186
187     //=================== Alice BODY parameters =============================
188     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
189
190
191     if (iMAG)
192     {
193         //=================== MAG parameters ============================
194         // --- Start with Magnet since detector layouts may be depending ---
195         // --- on the selected Magnet dimensions ---
196         AliMAG *MAG = new AliMAG("MAG", "Magnet");
197     }
198
199
200     if (iABSO)
201     {
202         //=================== ABSO parameters ============================
203         AliABSO *ABSO = new AliABSOv3("ABSO", "Muon Absorber");
204     }
205
206     if (iDIPO)
207     {
208         //=================== DIPO parameters ============================
209
210         AliDIPO *DIPO = new AliDIPOv3("DIPO", "Dipole version 3");
211     }
212
213     if (iHALL)
214     {
215         //=================== HALL parameters ============================
216
217         AliHALL *HALL = new AliHALLv3("HALL", "Alice Hall");
218     }
219
220
221     if (iFRAME)
222     {
223         //=================== FRAME parameters ============================
224
225         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
226         FRAME->SetHoles(1);
227     }
228
229     if (iSHIL)
230     {
231         //=================== SHIL parameters ============================
232
233         AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
234     }
235
236
237     if (iPIPE)
238     {
239         //=================== PIPE parameters ============================
240
241         AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
242     }
243  
244     if (iITS)
245     {
246         //=================== ITS parameters ============================
247
248         AliITS *ITS  = new AliITSv11Hybrid("ITS","ITS v11Hybrid");
249     }
250
251     if (iTPC)
252     {
253       //============================ TPC parameters =====================
254
255         AliTPC *TPC = new AliTPCv2("TPC", "Default");
256         TPC->SetPrimaryIonisation(); // not used with Geant3
257     }
258
259
260     if (iTOF) {
261         //=================== TOF parameters ============================
262
263         AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
264     }
265
266
267     if (iHMPID)
268     {
269         //=================== HMPID parameters ===========================
270
271         AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
272
273     }
274
275
276     if (iZDC)
277     {
278         //=================== ZDC parameters ============================
279
280         AliZDC *ZDC = new AliZDCv3("ZDC", "normal ZDC");
281     }
282
283     if (iTRD)
284     {
285         //=================== TRD parameters ============================
286
287         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
288         AliTRDgeometry *geoTRD = TRD->GetGeometry();
289         // Partial geometry: modules at 0,1,7,8,9,16,17
290         // starting at 3h in positive direction
291         geoTRD->SetSMstatus(2,0);
292         geoTRD->SetSMstatus(3,0);
293         geoTRD->SetSMstatus(4,0);
294         geoTRD->SetSMstatus(5,0);
295         geoTRD->SetSMstatus(6,0);
296         geoTRD->SetSMstatus(11,0);
297         geoTRD->SetSMstatus(12,0);
298         geoTRD->SetSMstatus(13,0);
299         geoTRD->SetSMstatus(14,0);
300         geoTRD->SetSMstatus(15,0);
301         geoTRD->SetSMstatus(16,0);
302     }
303
304     if (iFMD)
305     {
306         //=================== FMD parameters ============================
307
308         AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
309    }
310
311     if (iMUON)
312     {
313         //=================== MUON parameters ===========================
314         // New MUONv1 version (geometry defined via builders)
315
316         AliMUON *MUON = new AliMUONv1("MUON", "default");
317     }
318
319     if (iPHOS)
320     {
321         //=================== PHOS parameters ===========================
322
323        AliPHOS *PHOS = new AliPHOSv1("PHOS", "noCPV_Modules123");
324
325     }
326
327     if (iPMD)
328     {
329         //=================== PMD parameters ============================
330
331         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
332     }
333
334     if (iT0)
335     {
336         //=================== T0 parameters ============================
337         AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
338     }
339
340     if (iEMCAL)
341     {
342         //=================== EMCAL parameters ============================
343
344         AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_FIRSTYEAR");
345     }
346
347      if (iACORDE)
348     {
349         //=================== ACORDE parameters ============================
350
351         AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
352     }
353
354      if (iVZERO)
355     {
356         //=================== ACORDE parameters ============================
357
358         AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
359     }
360
361
362
363   // Load Geant4 + Geant4 VMC libraries
364   //
365   if (gClassTable->GetID("TGeant4") == -1) {
366     TString g4libsMacro = "$G4INSTALL/macro/g4libs.C";
367     //TString g4libsMacro = "$ALICE/geant4_vmc/examples/macro/g4libs.C";
368     // Load Geant4 libraries 
369     if (!gInterpreter->IsLoaded(g4libsMacro.Data())) {
370       gROOT->LoadMacro(g4libsMacro.Data());
371       gInterpreter->ProcessLine("g4libs()");
372     }
373   }    
374
375
376   // Create Geant4 VMC
377   //  
378   TGeant4 *geant4 = 0;
379   if ( ! gMC ) {
380     TG4RunConfiguration* runConfiguration=0x0;
381     for (Int_t iList = 0; iList < kListMax; iList++) {
382       if(iList<kListMax/2){
383         if(physicslist == iList){
384           runConfiguration = 
385             new TG4RunConfiguration("geomRoot", 
386                                     physicsListName[iList], 
387                                     "specialCuts+stackPopper+stepLimiter",
388                                     true);
389         }
390       }
391       else if(iList>=kListMax/2){//add "optical" PL to HadronPhysicsList
392         if(physicslist == iList){
393           runConfiguration = 
394             new TG4RunConfiguration("geomRoot", 
395                                     Form("%s+optical",physicsListName[iList-kListMax/2]), 
396                                     "specialCuts+stackPopper+stepLimiter",
397                                     true);
398         }
399       }
400     }
401     geant4 = new TGeant4("TGeant4", "The Geant4 Monte Carlo", runConfiguration);
402     cout << "Geant4 has been created." << endl;
403   } 
404   else {
405     cout << "Monte Carlo has been already created." << endl;
406   }  
407   
408   
409   
410   // Customization of Geant4 VMC
411   //
412   geant4->ProcessGeantCommand("/mcVerbose/all 1");  
413   geant4->ProcessGeantCommand("/mcVerbose/geometryManager 1");  
414   geant4->ProcessGeantCommand("/mcVerbose/opGeometryManager 1");  
415   geant4->ProcessGeantCommand("/mcTracking/loopVerbose 0");     
416   geant4->ProcessGeantCommand("/mcPhysics/rangeCuts 0.01 mm"); 
417   geant4->ProcessGeantCommand("/mcPhysics/selectOpProcess Scintillation");
418   geant4->ProcessGeantCommand("/mcPhysics/setOpProcessActivation false");
419   geant4->ProcessGeantCommand("/mcVerbose/composedPhysicsList 2");  
420   geant4->ProcessGeantCommand("/mcTracking/skipNeutrino true");
421   geant4->ProcessGeantCommand("/mcDet/setMaxStepInLowDensityMaterials 1 cm");
422
423   /*
424   // Set PAI model for TPC (TPC_Ne-CO2-N-2)
425   geant4->ProcessGeantCommand("/mcPhysics/emModel/selectMedium 219");
426   geant4->ProcessGeantCommand("/mcPhysics/emModel/setElossModel PAI");
427   geant4->ProcessGeantCommand("/mcPhysics/emModel/setFluctModel PAI");
428   geant4->ProcessGeantCommand("/mcPhysics/emModel/setParticles  all");
429   
430   // Set PAI model for TRD (TRD_XeCO2)
431   geant4->ProcessGeantCommand("/mcPhysics/emModel/selectMedium 291");
432   geant4->ProcessGeantCommand("/mcPhysics/emModel/setElossModel PAI");
433   geant4->ProcessGeantCommand("/mcPhysics/emModel/setFluctModel PAI");
434   geant4->ProcessGeantCommand("/mcPhysics/emModel/setParticles  all");
435   */
436
437
438  //
439   //=======================================================================
440   // ************* STEERING parameters FOR ALICE SIMULATION **************
441   // --- Specify event type to be tracked through the ALICE setup
442   // --- All positions are in cm, angles in degrees, and P and E in GeV
443
444
445     gMC->SetProcess("DCAY",1);
446     gMC->SetProcess("PAIR",1);
447     gMC->SetProcess("COMP",1);
448     gMC->SetProcess("PHOT",1);
449     gMC->SetProcess("PFIS",0);
450     gMC->SetProcess("DRAY",0);
451     gMC->SetProcess("ANNI",1);
452     gMC->SetProcess("BREM",1);
453     gMC->SetProcess("MUNU",1);
454     gMC->SetProcess("CKOV",1);
455     gMC->SetProcess("HADR",1);
456     gMC->SetProcess("LOSS",2);
457     gMC->SetProcess("MULS",1);
458     gMC->SetProcess("RAYL",1);
459
460     Float_t cut = 1.e-3;        // 1MeV cut by default
461     Float_t tofmax = 1.e10;
462
463     gMC->SetCut("CUTGAM", cut);
464     gMC->SetCut("CUTELE", cut);
465     gMC->SetCut("CUTNEU", cut);
466     gMC->SetCut("CUTHAD", cut);
467     gMC->SetCut("CUTMUO", cut);
468     gMC->SetCut("BCUTE",  cut); 
469     gMC->SetCut("BCUTM",  cut); 
470     gMC->SetCut("DCUTE",  cut); 
471     gMC->SetCut("DCUTM",  cut); 
472     gMC->SetCut("PPCUTM", cut);
473     gMC->SetCut("TOFMAX", tofmax); 
474
475
476
477
478   //======================//
479   // Set External decayer //
480   //======================//
481   TVirtualMCDecayer* decayer = new AliDecayerPythia();
482   decayer->SetForceDecay(kAll);
483   decayer->Init();
484   gMC->SetExternalDecayer(decayer);
485
486   //=========================//
487   // Generator Configuration //
488   //=========================//
489   AliGenerator* gener = 0x0;
490   
491   if (proc == kPythia6) {
492       gener = MbPythia();
493   } else if (proc == kPythia6D6T) {
494       gener = MbPythiaTuneD6T();
495   } else if (proc == kPythia6D6T_Flat) {
496       gener = MbPythiaTuneD6T_Flat();
497   } else if (proc == kPythia6ATLAS) {
498       gener = MbPythiaTuneATLAS();
499   } else if (proc == kPythia6ATLAS_Flat) {
500       gener = MbPythiaTuneATLAS_Flat();
501   } else if (proc == kPhojet) {
502       gener = MbPhojet();
503   } else if (proc == kPythia6_Perugia0) {
504       gener = MbPythiaTunePerugia0();
505   }
506
507   //
508   // PRIMARY VERTEX
509   //
510    Double_t xv = 0;
511    Double_t yv = 0;
512    Double_t zv = 0;
513
514    Double_t vtxPos[3];
515    Double_t vtxErr[3];
516    for (Int_t i = 0; i < 3; i++) {
517        vtxPos[i] = 0.;
518        vtxErr[i] = 0.;
519    }
520
521    if(AliCDBManager::Instance()->IsDefaultStorageSet()){
522       AliCDBEntry* entry = AliCDBManager::Instance()->Get("GRP/Calib/MeanVertexSPD");
523       AliESDVertex* vertex = dynamic_cast<AliESDVertex*> (entry->GetObject());
524       vertex->GetXYZ(vtxPos);
525       vertex->GetSigmaXYZ(vtxErr);
526     }
527
528     printf("Vertex position from OCDB entry: x = %13.3f, y = %13.3f, z = %13.3f (sigma = %13.3f)\n",
529     vtxPos[0], vtxPos[1], vtxPos[2], vtxErr[2]);
530
531      gener->SetOrigin(vtxPos[0], vtxPos[1], vtxPos[2]);   // vertex position
532   //
533   //
534   // Size of the interaction diamond
535   // Longitudinal
536      sigmaz = vtxErr[2];
537 //
538   // Transverse
539   Float_t betast  = 10;                 // beta* [m]
540   Float_t eps     = 2.5e-6;            // emittance [m]
541   Float_t gamma   = energy / 2.0 / 0.938272;  // relativistic gamma [1]
542   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.;  // [cm]
543  
544   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
545
546   gener->SetSigma(sigmaxy, sigmaxy, sigmaz);      // Sigma in (X,Y,Z) (cm) on IP position
547   gener->SetVertexSmear(kPerEvent);
548   gener->Init();
549
550
551
552
553
554
555
556
557
558
559
560
561
562 }
563 //
564 //           PYTHIA
565 //
566
567 AliGenerator* MbPythia()
568 {
569       comment = comment.Append(" pp: Pythia low-pt");
570 //
571 //    Pythia
572       AliGenPythia* pythia = new AliGenPythia(-1); 
573       pythia->SetMomentumRange(0, 999999.);
574       pythia->SetThetaRange(0., 180.);
575       pythia->SetYRange(-12.,12.);
576       pythia->SetPtRange(0,1000.);
577       pythia->SetProcess(kPyMb);
578       pythia->SetEnergyCMS(energy);
579       
580       return pythia;
581 }
582
583 AliGenerator* MbPythiaTuneD6T()
584 {
585       comment = comment.Append(" pp: Pythia low-pt");
586 //
587 //    Pythia
588       AliGenPythia* pythia = new AliGenPythia(-1); 
589       pythia->SetMomentumRange(0, 999999.);
590       pythia->SetThetaRange(0., 180.);
591       pythia->SetYRange(-12.,12.);
592       pythia->SetPtRange(0,1000.);
593       pythia->SetProcess(kPyMb);
594       pythia->SetEnergyCMS(energy);
595 //    Tune
596 //    109     D6T : Rick Field's CDF Tune D6T (NB: needs CTEQ6L pdfs externally)
597       pythia->SetTune(109); // F I X 
598       pythia->SetStrucFunc(kCTEQ6l);
599 //
600       return pythia;
601 }
602
603 AliGenerator* MbPythiaTunePerugia0()
604 {
605       comment = comment.Append(" pp: Pythia Minimum Bias Tune Perugia0");
606 //
607 //    Pythia
608       AliGenPythia* pythia = new AliGenPythia(-1); 
609       pythia->SetMomentumRange(0, 999999.);
610       pythia->SetThetaRange(0., 180.);
611       pythia->SetYRange(-12.,12.);
612       pythia->SetPtRange(0, 1000.);
613       pythia->SetProcess(kPyMb);
614       pythia->SetEnergyCMS(energy);
615 //    Tune
616       pythia->SetTune(320);
617 //
618       return pythia;
619 }
620
621 AliGenerator* MbPythiaTuneATLAS()
622 {
623       comment = comment.Append(" pp: Pythia low-pt");
624 //
625 //    Pythia
626       AliGenPythia* pythia = new AliGenPythia(-1); 
627       pythia->SetMomentumRange(0, 999999.);
628       pythia->SetThetaRange(0., 180.);
629       pythia->SetYRange(-12.,12.);
630       pythia->SetPtRange(0,1000.);
631       pythia->SetProcess(kPyMb);
632       pythia->SetEnergyCMS(energy);
633 //    Tune
634 //    C   306 ATLAS-CSC: Arthur Moraes' (new) ATLAS tune (needs CTEQ6L externally)
635       pythia->SetTune(306);
636       pythia->SetStrucFunc(kCTEQ6l);
637 //
638       return pythia;
639 }
640
641 AliGenerator* MbPythiaTuneATLAS_Flat()
642 {
643       AliGenPythia* pythia = MbPythiaTuneATLAS();
644       
645       comment = comment.Append("; flat multiplicity distribution");
646       
647       // set high multiplicity trigger
648       // this weight achieves a flat multiplicity distribution
649       TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
650       weight->SetBinContent(1,5.49443);
651       weight->SetBinContent(2,8.770816);
652       weight->SetBinContent(6,0.4568624);
653       weight->SetBinContent(7,0.2919915);
654       weight->SetBinContent(8,0.6674189);
655       weight->SetBinContent(9,0.364737);
656       weight->SetBinContent(10,0.8818444);
657       weight->SetBinContent(11,0.531885);
658       weight->SetBinContent(12,1.035197);
659       weight->SetBinContent(13,0.9394057);
660       weight->SetBinContent(14,0.9643193);
661       weight->SetBinContent(15,0.94543);
662       weight->SetBinContent(16,0.9426507);
663       weight->SetBinContent(17,0.9423649);
664       weight->SetBinContent(18,0.789456);
665       weight->SetBinContent(19,1.149026);
666       weight->SetBinContent(20,1.100491);
667       weight->SetBinContent(21,0.6350525);
668       weight->SetBinContent(22,1.351941);
669       weight->SetBinContent(23,0.03233504);
670       weight->SetBinContent(24,0.9574557);
671       weight->SetBinContent(25,0.868133);
672       weight->SetBinContent(26,1.030998);
673       weight->SetBinContent(27,1.08897);
674       weight->SetBinContent(28,1.251382);
675       weight->SetBinContent(29,0.1391099);
676       weight->SetBinContent(30,1.192876);
677       weight->SetBinContent(31,0.448944);
678       weight->SetBinContent(32,1);
679       weight->SetBinContent(33,1);
680       weight->SetBinContent(34,1);
681       weight->SetBinContent(35,1);
682       weight->SetBinContent(36,0.9999997);
683       weight->SetBinContent(37,0.9999997);
684       weight->SetBinContent(38,0.9999996);
685       weight->SetBinContent(39,0.9999996);
686       weight->SetBinContent(40,0.9999995);
687       weight->SetBinContent(41,0.9999993);
688       weight->SetBinContent(42,1);
689       weight->SetBinContent(43,1);
690       weight->SetBinContent(44,1);
691       weight->SetBinContent(45,1);
692       weight->SetBinContent(46,1);
693       weight->SetBinContent(47,0.9999999);
694       weight->SetBinContent(48,0.9999998);
695       weight->SetBinContent(49,0.9999998);
696       weight->SetBinContent(50,0.9999999);
697       weight->SetBinContent(51,0.9999999);
698       weight->SetBinContent(52,0.9999999);
699       weight->SetBinContent(53,0.9999999);
700       weight->SetBinContent(54,0.9999998);
701       weight->SetBinContent(55,0.9999998);
702       weight->SetBinContent(56,0.9999998);
703       weight->SetBinContent(57,0.9999997);
704       weight->SetBinContent(58,0.9999996);
705       weight->SetBinContent(59,0.9999995);
706       weight->SetBinContent(60,1);
707       weight->SetBinContent(61,1);
708       weight->SetBinContent(62,1);
709       weight->SetBinContent(63,1);
710       weight->SetBinContent(64,1);
711       weight->SetBinContent(65,0.9999999);
712       weight->SetBinContent(66,0.9999998);
713       weight->SetBinContent(67,0.9999998);
714       weight->SetBinContent(68,0.9999999);
715       weight->SetBinContent(69,1);
716       weight->SetBinContent(70,1);
717       weight->SetBinContent(71,0.9999997);
718       weight->SetBinContent(72,0.9999995);
719       weight->SetBinContent(73,0.9999994);
720       weight->SetBinContent(74,1);
721       weight->SetBinContent(75,1);
722       weight->SetBinContent(76,1);
723       weight->SetBinContent(77,1);
724       weight->SetBinContent(78,0.9999999);
725       weight->SetBinContent(79,1);
726       weight->SetBinContent(80,1);
727       weight->SetEntries(526);
728         
729       Int_t limit = weight->GetRandom();
730       pythia->SetTriggerChargedMultiplicity(limit, 1.4);
731       
732       comment = comment.Append(Form("; multiplicity threshold set to %d in |eta| < 1.4", limit));
733
734       return pythia;
735 }
736 AliGenerator* MbPythiaTuneD6T_Flat()
737 {
738       AliGenPythia* pythia = MbPythiaTuneD6T();
739       
740       comment = comment.Append("; flat multiplicity distribution");
741       
742       // set high multiplicity trigger
743       // this weight achieves a flat multiplicity distribution
744       TH1 *weight = new TH1D("weight","weight",201,-0.5,200.5);
745       weight->SetBinContent(1,5.49443);
746       weight->SetBinContent(2,8.770816);
747       weight->SetBinContent(6,0.4568624);
748       weight->SetBinContent(7,0.2919915);
749       weight->SetBinContent(8,0.6674189);
750       weight->SetBinContent(9,0.364737);
751       weight->SetBinContent(10,0.8818444);
752       weight->SetBinContent(11,0.531885);
753       weight->SetBinContent(12,1.035197);
754       weight->SetBinContent(13,0.9394057);
755       weight->SetBinContent(14,0.9643193);
756       weight->SetBinContent(15,0.94543);
757       weight->SetBinContent(16,0.9426507);
758       weight->SetBinContent(17,0.9423649);
759       weight->SetBinContent(18,0.789456);
760       weight->SetBinContent(19,1.149026);
761       weight->SetBinContent(20,1.100491);
762       weight->SetBinContent(21,0.6350525);
763       weight->SetBinContent(22,1.351941);
764       weight->SetBinContent(23,0.03233504);
765       weight->SetBinContent(24,0.9574557);
766       weight->SetBinContent(25,0.868133);
767       weight->SetBinContent(26,1.030998);
768       weight->SetBinContent(27,1.08897);
769       weight->SetBinContent(28,1.251382);
770       weight->SetBinContent(29,0.1391099);
771       weight->SetBinContent(30,1.192876);
772       weight->SetBinContent(31,0.448944);
773       weight->SetBinContent(32,1);
774       weight->SetBinContent(33,1);
775       weight->SetBinContent(34,1);
776       weight->SetBinContent(35,1);
777       weight->SetBinContent(36,0.9999997);
778       weight->SetBinContent(37,0.9999997);
779       weight->SetBinContent(38,0.9999996);
780       weight->SetBinContent(39,0.9999996);
781       weight->SetBinContent(40,0.9999995);
782       weight->SetBinContent(41,0.9999993);
783       weight->SetBinContent(42,1);
784       weight->SetBinContent(43,1);
785       weight->SetBinContent(44,1);
786       weight->SetBinContent(45,1);
787       weight->SetBinContent(46,1);
788       weight->SetBinContent(47,0.9999999);
789       weight->SetBinContent(48,0.9999998);
790       weight->SetBinContent(49,0.9999998);
791       weight->SetBinContent(50,0.9999999);
792       weight->SetBinContent(51,0.9999999);
793       weight->SetBinContent(52,0.9999999);
794       weight->SetBinContent(53,0.9999999);
795       weight->SetBinContent(54,0.9999998);
796       weight->SetBinContent(55,0.9999998);
797       weight->SetBinContent(56,0.9999998);
798       weight->SetBinContent(57,0.9999997);
799       weight->SetBinContent(58,0.9999996);
800       weight->SetBinContent(59,0.9999995);
801       weight->SetBinContent(60,1);
802       weight->SetBinContent(61,1);
803       weight->SetBinContent(62,1);
804       weight->SetBinContent(63,1);
805       weight->SetBinContent(64,1);
806       weight->SetBinContent(65,0.9999999);
807       weight->SetBinContent(66,0.9999998);
808       weight->SetBinContent(67,0.9999998);
809       weight->SetBinContent(68,0.9999999);
810       weight->SetBinContent(69,1);
811       weight->SetBinContent(70,1);
812       weight->SetBinContent(71,0.9999997);
813       weight->SetBinContent(72,0.9999995);
814       weight->SetBinContent(73,0.9999994);
815       weight->SetBinContent(74,1);
816       weight->SetBinContent(75,1);
817       weight->SetBinContent(76,1);
818       weight->SetBinContent(77,1);
819       weight->SetBinContent(78,0.9999999);
820       weight->SetBinContent(79,1);
821       weight->SetBinContent(80,1);
822       weight->SetEntries(526);
823         
824       Int_t limit = weight->GetRandom();
825       pythia->SetTriggerChargedMultiplicity(limit, 1.4);
826       
827       comment = comment.Append(Form("; multiplicity threshold set to %d in |eta| < 1.4", limit));
828
829       return pythia;
830 }
831
832 AliGenerator* MbPhojet()
833 {
834       comment = comment.Append(" pp: Pythia low-pt");
835 //
836 //    DPMJET
837 #if defined(__CINT__)
838   gSystem->Load("libdpmjet");      // Parton density functions
839   gSystem->Load("libTDPMjet");      // Parton density functions
840 #endif
841       AliGenDPMjet* dpmjet = new AliGenDPMjet(-1); 
842       dpmjet->SetMomentumRange(0, 999999.);
843       dpmjet->SetThetaRange(0., 180.);
844       dpmjet->SetYRange(-12.,12.);
845       dpmjet->SetPtRange(0,1000.);
846       dpmjet->SetProcess(kDpmMb);
847       dpmjet->SetEnergyCMS(energy);
848
849       return dpmjet;
850 }
851
852 void ProcessEnvironmentVars()
853 {
854     // Run type
855     if (gSystem->Getenv("CONFIG_RUN_TYPE")) {
856       for (Int_t iRun = 0; iRun < kRunMax; iRun++) {
857         if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun])==0) {
858           proc = (PDC06Proc_t)iRun;
859           cout<<"Run type set to "<<pprRunName[iRun]<<endl;
860         }
861       }
862     }
863
864     // Field
865     if (gSystem->Getenv("CONFIG_FIELD")) {
866       for (Int_t iField = 0; iField < kFieldMax; iField++) {
867         if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField])==0) {
868           mag = (Mag_t)iField;
869           cout<<"Field set to "<<pprField[iField]<<endl;
870         }
871       }
872     }
873
874     // Energy
875     if (gSystem->Getenv("CONFIG_ENERGY")) {
876       energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
877       cout<<"Energy set to "<<energy<<" GeV"<<endl;
878     }
879
880     // Random Number seed
881     if (gSystem->Getenv("CONFIG_SEED")) {
882       seed = atoi(gSystem->Getenv("CONFIG_SEED"));
883     }
884
885     // Physics lists
886     if (gSystem->Getenv("CONFIG_PHYSICSLIST")) {
887       for (Int_t iList = 0; iList < kListMax; iList++) {
888         if (strcmp(gSystem->Getenv("CONFIG_PHYSICSLIST"), physicsListName[iList])==0) {
889           physicslist = (PhysicsList_t)iList;
890           cout<<"Physics list set to "<< physicsListName[iList]<<endl;
891         }
892       }
893     }
894
895 }