update for new geant4 production
[u/mrichter/AliRoot.git] / test / vmctest / production / Config.C
CommitLineData
ba2958c4 1// $Id$
b2b18d9e 2//
ba2958c4 3// Configuration for the Geant4 production 2010
4// By E. Sicking, CERN
b2b18d9e 5//
6
b2b18d9e 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
48enum PDC06Proc_t
49{
50 kPythia6, kPythia6D6T, kPythia6D6T_Flat, kPythia6ATLAS, kPythia6ATLAS_Flat, kPythia6_Perugia0, kPhojet, kRunMax
51};
52
53const char * pprRunName[] = {
54 "kPythia6", "kPythia6D6T", "kPythia6D6T_Flat", "kPythia6ATLAS", "kPythia6ATLAS_Flat", "kPythia6_Perugia0", "kPhojet"
55};
56
57enum Mag_t
58{
59 kNoField, k5kG, kFieldMax
60};
61
62const char * pprField[] = {
63 "kNoField", "k5kG"
64};
65
66enum PhysicsList_t
67{
68 QGSP_BERT_EMV,CHIPS,QGSP_BERT_CHIPS,QGSP_BERT_EMV_OPTICAL, CHIPS_OPTICAL, QGSP_BERT_CHIPS_OPTICAL, kListMax
69};
70
71const char * physicsListName[] = {
72 "QGSP_BERT_EMV", "CHIPS", "QGSP_BERT_CHIPS",
73 "QGSP_BERT_EMV_OPTICAL", "CHIPS_OPTICAL", "QGSP_BERT_CHIPS_OPTICAL"
74};
75
6b4f7fc0 76enum PprTrigConf_t
77{
78 kDefaultPPTrig, kDefaultPbPbTrig
79};
80
81const char * pprTrigConfName[] = {
82 "p-p","Pb-Pb"
83};
84
b2b18d9e 85//--- Functions ---
86class AliGenPythia;
87AliGenerator *MbPythia();
88AliGenerator *MbPythiaTuneD6T();
89AliGenerator *MbPhojet();
90void ProcessEnvironmentVars();
91
92// Geterator, field, beam energy
93static PDC06Proc_t proc = kPhojet;
94static Mag_t mag = k5kG;
95static Float_t energy = 10000; // energy in CMS
96static PhysicsList_t physicslist = QGSP_BERT_EMV;
6b4f7fc0 97static PprTrigConf_t strig = kDefaultPPTrig; // default pp trigger configuration
98
b2b18d9e 99//========================//
100// Set Random Number seed //
101//========================//
102TDatime dt;
103static UInt_t seed = dt.Get();
104
105// Comment line
106static TString comment;
107
108void 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
6b4f7fc0 157 AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
158 cout<<"Trigger configuration is set to "<<pprTrigConfName[strig]<<endl;
b2b18d9e 159
b2b18d9e 160 printf("\n \n Comment: %s \n \n", comment.Data());
fd03ab3d 161
b2b18d9e 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");
c9ec33c4 256 TPC->SetPrimaryIonisation(); // not used with Geant3
b2b18d9e 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
fd03ab3d 323 AliPHOS *PHOS = new AliPHOSv1("PHOS", "noCPV_Modules123");
b2b18d9e 324
b2b18d9e 325 }
326
b2b18d9e 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) {
ba2958c4 366 TString g4libsMacro = "$G4INSTALL/macro/g4libs.C";
367 //TString g4libsMacro = "$ALICE/geant4_vmc/examples/macro/g4libs.C";
b2b18d9e 368 // Load Geant4 libraries
ba2958c4 369 if (!gInterpreter->IsLoaded(g4libsMacro.Data())) {
370 gROOT->LoadMacro(g4libsMacro.Data());
b2b18d9e 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");
e6918581 420 geant4->ProcessGeantCommand("/mcTracking/skipNeutrino true");
c9ec33c4 421 geant4->ProcessGeantCommand("/mcDet/setMaxStepInLowDensityMaterials 1 cm");
b2b18d9e 422
fd03ab3d 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 */
b2b18d9e 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
567AliGenerator* 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
583AliGenerator* 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
603AliGenerator* 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
621AliGenerator* 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
641AliGenerator* 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}
736AliGenerator* 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
832AliGenerator* 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
852void 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}