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