test/generators/herwig/Config.C

   1 //
   2 // Configuration for the first physics production 2008
   3 //
   4
   5 // One can use the configuration macro in compiled mode by
   6 // root [0] gSystem->Load("libgeant321");
   7 // root [0] gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/include\
   8 //                   -I$ALICE_ROOT -I$ALICE/geant3/TGeant3");
   9 // root [0] .x grun.C(1,"Config.C++")
  10
  11 #if !defined(__CINT__) || defined(__MAKECINT__)
  12 #include <Riostream.h>
  13 #include <TRandom.h>
  14 #include <TDatime.h>
  15 #include <TSystem.h>
  16 #include <TVirtualMC.h>
  17 #include <TGeant3TGeo.h>
  18 #include "STEER/AliRunLoader.h"
  19 #include "STEER/AliRun.h"
  20 #include "STEER/AliConfig.h"
  21 #include "PYTHIA6/AliDecayerPythia.h"
  22 #include "PYTHIA6/AliGenPythia.h"
  23 #include "TDPMjet/AliGenDPMjet.h"
  24 #include "THerwig/AliGenHerwig.h"
  25 #include "STEER/AliMagFCheb.h"
  26 #include "STRUCT/AliBODY.h"
  27 #include "STRUCT/AliMAG.h"
  28 #include "STRUCT/AliABSOv3.h"
  29 #include "STRUCT/AliDIPOv3.h"
  30 #include "STRUCT/AliHALLv3.h"
  31 #include "STRUCT/AliFRAMEv2.h"
  32 #include "STRUCT/AliSHILv3.h"
  33 #include "STRUCT/AliPIPEv3.h"
  34 #include "ITS/AliITSv11.h"
  35 #include "TPC/AliTPCv2.h"
  36 #include "TOF/AliTOFv6T0.h"
  37 #include "HMPID/AliHMPIDv3.h"
  38 #include "ZDC/AliZDCv4.h"
  39 #include "TRD/AliTRDv1.h"
  40 #include "TRD/AliTRDgeometry.h"
  41 #include "FMD/AliFMDv1.h"
  42 #include "MUON/AliMUONv1.h"
  43 #include "PHOS/AliPHOSv1.h"
  44 #include "PHOS/AliPHOSSimParam.h"
  45 #include "PMD/AliPMDv1.h"
  46 #include "T0/AliT0v1.h"
  47 #include "EMCAL/AliEMCALv2.h"
  48 #include "ACORDE/AliACORDEv1.h"
  49 #include "VZERO/AliVZEROv7.h"
  50 #endif
  51
  52 enum PDC06Proc_t
  53 {
  54   kPythia6, kPhojet, kHerwig, kRunMax
  55 };
  56
  57 const char * pprRunName[] =
  58 { "kPythia6", "kPhojet", "kHerwig" };
  59
  60 enum Mag_t
  61 {
  62   kNoField, k5kG, kFieldMax
  63 };
  64
  65 const char * pprField[] =
  66 { "kNoField", "k5kG" };
  67
  68 enum PprTrigConf_t
  69 {
  70   kDefaultPPTrig, kDefaultPbPbTrig
  71 };
  72
  73 const char * pprTrigConfName[] =
  74 { "p-p", "Pb-Pb" };
  75
  76 static PprTrigConf_t strig = kDefaultPPTrig;// default PP trigger configuration
  77
  78 //--- Functions ---
  79
  80 class AliGenPythia;
  81 AliGenerator *MbPythia();
  82 AliGenerator *MbPhojet();
  83 AliGenerator *Herwig();
  84 void ProcessEnvironmentVars();
  85
  86 // Geterator, field, beam energy
  87 static PDC06Proc_t proc = kHerwig;
  88 static Mag_t mag = k5kG;
  89 static Float_t energy = 10000; // energy in CMS
  90 //========================//
  91 // Set Random Number seed //
  92 //========================//
  93 TDatime dt;
  94 static UInt_t seed = dt.Get();
  95
  96 // Comment line
  97 static TString comment;
  98
  99 void Config()
 100 {
 101
 102   // Get settings from environment variables
 103   ProcessEnvironmentVars();
 104
 105   gRandom->SetSeed(seed);
 106   cerr << "Seed for random number generation= " << seed << endl;
 107
 108   // Libraries required by geant321
 109 #if defined(__CINT__)
 110   gSystem->Load("liblhapdf"); // Parton density functions
 111   gSystem->Load("libEGPythia6"); // TGenerator interface
 112   gSystem->Load("libpythia6"); // Pythia
 113   gSystem->Load("libAliPythia6"); // ALICE specific implementations
 114   gSystem->Load("libgeant321");
 115 #endif
 116
 117   new TGeant3TGeo("C++ Interface to Geant3");
 118
 119   //=======================================================================
 120   //  Create the output file
 121
 122
 123   AliRunLoader* rl = 0x0;
 124
 125   cout << "Config.C: Creating Run Loader ..." << endl;
 126   rl = AliRunLoader::Open("galice.root", AliConfig::GetDefaultEventFolderName(), "recreate");
 127   if (rl == 0x0)
 128   {
 129     gAlice->Fatal("Config.C", "Can not instatiate the Run Loader");
 130     return;
 131   }
 132   rl->SetCompressionLevel(2);
 133   rl->SetNumberOfEventsPerFile(1000);
 134   gAlice->SetRunLoader(rl);
 135
 136   // Set the trigger configuration: proton-proton
 137   AliSimulation::Instance()->SetTriggerConfig(pprTrigConfName[strig]);
 138   cout << "Trigger configuration is set to  " << pprTrigConfName[strig] << endl;
 139
 140   //
 141   //=======================================================================
 142   // ************* STEERING parameters FOR ALICE SIMULATION **************
 143   // --- Specify event type to be tracked through the ALICE setup
 144   // --- All positions are in cm, angles in degrees, and P and E in GeV
 145
 146
 147   gMC->SetProcess("DCAY", 1);
 148   gMC->SetProcess("PAIR", 1);
 149   gMC->SetProcess("COMP", 1);
 150   gMC->SetProcess("PHOT", 1);
 151   gMC->SetProcess("PFIS", 0);
 152   gMC->SetProcess("DRAY", 0);
 153   gMC->SetProcess("ANNI", 1);
 154   gMC->SetProcess("BREM", 1);
 155   gMC->SetProcess("MUNU", 1);
 156   gMC->SetProcess("CKOV", 1);
 157   gMC->SetProcess("HADR", 1);
 158   gMC->SetProcess("LOSS", 2);
 159   gMC->SetProcess("MULS", 1);
 160   gMC->SetProcess("RAYL", 1);
 161
 162   Float_t cut = 1.e-3; // 1MeV cut by default
 163   Float_t tofmax = 1.e10;
 164
 165   gMC->SetCut("CUTGAM", cut);
 166   gMC->SetCut("CUTELE", cut);
 167   gMC->SetCut("CUTNEU", cut);
 168   gMC->SetCut("CUTHAD", cut);
 169   gMC->SetCut("CUTMUO", cut);
 170   gMC->SetCut("BCUTE", cut);
 171   gMC->SetCut("BCUTM", cut);
 172   gMC->SetCut("DCUTE", cut);
 173   gMC->SetCut("DCUTM", cut);
 174   gMC->SetCut("PPCUTM", cut);
 175   gMC->SetCut("TOFMAX", tofmax);
 176
 177   //======================//
 178   // Set External decayer //
 179   //======================//
 180   TVirtualMCDecayer* decayer = new AliDecayerPythia();
 181   decayer->SetForceDecay(kAll);
 182   decayer->Init();
 183   gMC->SetExternalDecayer(decayer);
 184
 185   //=========================//
 186   // Generator Configuration //
 187   //=========================//
 188   AliGenerator* gener = 0x0;
 189
 190   switch (proc)
 191   {
 192   case kPythia6:
 193     gener = MbPythia();
 194     break;
 195   case kPhojet:
 196     gener = MbPhojet();
 197     break;
 198   case kHerwig:
 199     gener = Herwig();
 200     break;
 201   }
 202
 203   // PRIMARY VERTEX
 204   //
 205   gener->SetOrigin(0., 0., 0.); // vertex position
 206   //
 207   // Size of the interaction diamond
 208   // Longitudinal
 209   Float_t sigmaz = 5.4 / TMath::Sqrt(2.); // [cm]
 210   if (energy == 900)
 211     sigmaz = 10.5 / TMath::Sqrt(2.); // [cm]
 212   //
 213   // Transverse
 214   Float_t betast = 10; // beta* [m]
 215   Float_t eps = 3.75e-6; // emittance [m]
 216   Float_t gamma = energy / 2.0 / 0.938272; // relativistic gamma [1]
 217   Float_t sigmaxy = TMath::Sqrt(eps * betast / gamma) / TMath::Sqrt(2.) * 100.; // [cm]
 218   printf("\n \n Diamond size x-y: %10.3e z: %10.3e\n \n", sigmaxy, sigmaz);
 219
 220   gener->SetSigma(sigmaxy, sigmaxy, sigmaz); // Sigma in (X,Y,Z) (cm) on IP position
 221   gener->SetCutVertexZ(3.); // Truncate at 3 sigma
 222   gener->SetVertexSmear(kPerEvent);
 223
 224   gener->Init();
 225
 226   // FIELD
 227   AliMagF* field = 0x0;
 228
 229   if (mag == kNoField)
 230   {
 231     comment = comment.Append(" | L3 field 0.0 T");
 232     field = new AliMagF("Maps", "Maps", -1., -1., AliMagF::k2kG);
 233   }
 234   else if (mag == k5kG)
 235   {
 236     comment = comment.Append(" | L3 field 0.5 T");
 237     field = new AliMagF("Maps", "Maps", -1., -1., AliMagF::k5kG);
 238   }
 239   printf("\n \n Comment: %s \n \n", comment.Data());
 240
 241   TGeoGlobalMagField::Instance()->SetField(field);
 242
 243   rl->CdGAFile();
 244
 245   Int_t iABSO = 1;
 246   Int_t iACORDE = 0;
 247   Int_t iDIPO = 1;
 248   Int_t iEMCAL = 0;
 249   Int_t iFMD = 1;
 250   Int_t iFRAME = 1;
 251   Int_t iHALL = 1;
 252   Int_t iITS = 1;
 253   Int_t iMAG = 1;
 254   Int_t iMUON = 1;
 255   Int_t iPHOS = 1;
 256   Int_t iPIPE = 1;
 257   Int_t iPMD = 0;
 258   Int_t iHMPID = 1;
 259   Int_t iSHIL = 1;
 260   Int_t iT0 = 1;
 261   Int_t iTOF = 1;
 262   Int_t iTPC = 1;
 263   Int_t iTRD = 1;
 264   Int_t iVZERO = 1;
 265   Int_t iZDC = 1;
 266
 267   //=================== Alice BODY parameters =============================
 268   AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
 269
 270   if (iMAG)
 271   {
 272     //=================== MAG parameters ============================
 273     // --- Start with Magnet since detector layouts may be depending ---
 274     // --- on the selected Magnet dimensions ---
 275     AliMAG *MAG = new AliMAG("MAG", "Magnet");
 276   }
 277
 278   if (iABSO)
 279   {
 280     //=================== ABSO parameters ============================
 281     AliABSO *ABSO = new AliABSOv3("ABSO", "Muon Absorber");
 282   }
 283
 284   if (iDIPO)
 285   {
 286     //=================== DIPO parameters ============================
 287
 288     AliDIPO *DIPO = new AliDIPOv3("DIPO", "Dipole version 3");
 289   }
 290
 291   if (iHALL)
 292   {
 293     //=================== HALL parameters ============================
 294
 295     AliHALL *HALL = new AliHALLv3("HALL", "Alice Hall");
 296   }
 297
 298   if (iFRAME)
 299   {
 300     //=================== FRAME parameters ============================
 301
 302     AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
 303     FRAME->SetHoles(1);
 304   }
 305
 306   if (iSHIL)
 307   {
 308     //=================== SHIL parameters ============================
 309
 310     AliSHIL *SHIL = new AliSHILv3("SHIL", "Shielding Version 3");
 311   }
 312
 313   if (iPIPE)
 314   {
 315     //=================== PIPE parameters ============================
 316
 317     AliPIPE *PIPE = new AliPIPEv3("PIPE", "Beam Pipe");
 318   }
 319
 320   if (iITS)
 321   {
 322     //=================== ITS parameters ============================
 323
 324     AliITS *ITS = new AliITSv11("ITS", "ITS v11");
 325   }
 326
 327   if (iTPC)
 328   {
 329     //============================ TPC parameters =====================
 330
 331     AliTPC *TPC = new AliTPCv2("TPC", "Default");
 332   }
 333
 334   if (iTOF)
 335   {
 336     //=================== TOF parameters ============================
 337
 338     AliTOF *TOF = new AliTOFv6T0("TOF", "normal TOF");
 339   }
 340
 341   if (iHMPID)
 342   {
 343     //=================== HMPID parameters ===========================
 344
 345     AliHMPID *HMPID = new AliHMPIDv3("HMPID", "normal HMPID");
 346
 347   }
 348
 349   if (iZDC)
 350   {
 351     //=================== ZDC parameters ============================
 352
 353     AliZDC *ZDC = new AliZDCv4("ZDC", "normal ZDC");
 354   }
 355
 356   if (iTRD)
 357   {
 358     //=================== TRD parameters ============================
 359
 360     AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
 361     AliTRDgeometry *geoTRD = TRD->GetGeometry();
 362     // Partial geometry: modules at 0,1,7,8,9,10,17
 363     // starting at 3h in positive direction
 364     geoTRD->SetSMstatus(2, 0);
 365     geoTRD->SetSMstatus(3, 0);
 366     geoTRD->SetSMstatus(4, 0);
 367     geoTRD->SetSMstatus(5, 0);
 368     geoTRD->SetSMstatus(6, 0);
 369     geoTRD->SetSMstatus(11, 0);
 370     geoTRD->SetSMstatus(12, 0);
 371     geoTRD->SetSMstatus(13, 0);
 372     geoTRD->SetSMstatus(14, 0);
 373     geoTRD->SetSMstatus(15, 0);
 374     geoTRD->SetSMstatus(16, 0);
 375   }
 376
 377   if (iFMD)
 378   {
 379     //=================== FMD parameters ============================
 380
 381     AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
 382   }
 383
 384   if (iMUON)
 385   {
 386     //=================== MUON parameters ===========================
 387     // New MUONv1 version (geometry defined via builders)
 388
 389     AliMUON *MUON = new AliMUONv1("MUON", "default");
 390   }
 391
 392   if (iPHOS)
 393   {
 394     //=================== PHOS parameters ===========================
 395
 396     AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
 397   }
 398
 399   if (iPMD)
 400   {
 401     //=================== PMD parameters ============================
 402
 403     AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
 404   }
 405
 406   if (iT0)
 407   {
 408     //=================== T0 parameters ============================
 409     AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
 410   }
 411
 412   if (iEMCAL)
 413   {
 414     //=================== EMCAL parameters ============================
 415
 416     AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "EMCAL_COMPLETEV1");
 417   }
 418
 419   if (iACORDE)
 420   {
 421     //=================== ACORDE parameters ============================
 422
 423     AliACORDE *ACORDE = new AliACORDEv1("ACORDE", "normal ACORDE");
 424   }
 425
 426   if (iVZERO)
 427   {
 428     //=================== ACORDE parameters ============================
 429
 430     AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
 431   }
 432 }
 433
 434 //           PYTHIA
 435 AliGenerator* MbPythia()
 436 {
 437   comment = comment.Append(" pp at 14 TeV: Pythia low-pt");
 438
 439   AliGenPythia* pythia = new AliGenPythia(-1);
 440   pythia->SetMomentumRange(0, 999999.);
 441   pythia->SetThetaRange(0., 180.);
 442   pythia->SetYRange(-12., 12.);
 443   pythia->SetPtRange(0, 1000.);
 444   pythia->SetProcess(kPyMb);
 445   pythia->SetEnergyCMS(energy);
 446
 447   return pythia;
 448 }
 449
 450 //          DPMJET
 451 AliGenerator* MbPhojet()
 452 {
 453   comment = comment.Append(" pp at 14 TeV: Phojet low-pt");
 454
 455 #if defined(__CINT__)
 456   gSystem->Load("libdpmjet"); // Parton density functions
 457   gSystem->Load("libTDPMjet"); // Parton density functions
 458 #endif
 459
 460   AliGenDPMjet* dpmjet = new AliGenDPMjet(-1);
 461   dpmjet->SetMomentumRange(0, 999999.);
 462   dpmjet->SetThetaRange(0., 180.);
 463   dpmjet->SetYRange(-12., 12.);
 464   dpmjet->SetPtRange(0, 1000.);
 465   dpmjet->SetProcess(kDpmMb);
 466   dpmjet->SetEnergyCMS(energy);
 467
 468   return dpmjet;
 469 }
 470
 471 //          HERWIG
 472 AliGenerator* Herwig()
 473 {
 474   comment = comment.Append("pp at 14 TeV: Herwig");
 475
 476 #if defined(__CINT__)
 477   gSystem->Load("libherwig"); // HERWIG library
 478   gSystem->Load("libTHerwig"); // HERWIG library
 479 #endif
 480
 481     AliGenHerwig *herwig = new AliGenHerwig(-1);
 482     // final state kinematic cuts
 483     herwig->SetOrigin(0,0,0); // vertex position
 484     herwig->SetVertexSmear(kPerEvent);
 485     herwig->SetSigma(0, 0, 5.6); // Sigma in (X,Y,Z) (cm) on IP position
 486     // Beams
 487     herwig->SetProjectile("P");
 488     herwig->SetTarget("P");
 489     // Beam momenta
 490     herwig->SetBeamMomenta(3500., 3500.);
 491     // bbar
 492     herwig->SetProcess(1705);
 493
 494   return herwig;
 495 }
 496
 497 void ProcessEnvironmentVars()
 498 {
 499   // Run type
 500   if (gSystem->Getenv("CONFIG_RUN_TYPE"))
 501   {
 502     for (Int_t iRun = 0; iRun < kRunMax; iRun++)
 503     {
 504       if (strcmp(gSystem->Getenv("CONFIG_RUN_TYPE"), pprRunName[iRun]) == 0)
 505       {
 506         proc = (PDC06Proc_t) iRun;
 507         cout << "Run type set to " << pprRunName[iRun] << endl;
 508       }
 509     }
 510   }
 511
 512   // Field
 513   if (gSystem->Getenv("CONFIG_FIELD"))
 514   {
 515     for (Int_t iField = 0; iField < kFieldMax; iField++)
 516     {
 517       if (strcmp(gSystem->Getenv("CONFIG_FIELD"), pprField[iField]) == 0)
 518       {
 519         mag = (Mag_t) iField;
 520         cout << "Field set to " << pprField[iField] << endl;
 521       }
 522     }
 523   }
 524
 525   // Energy
 526   if (gSystem->Getenv("CONFIG_ENERGY"))
 527   {
 528     energy = atoi(gSystem->Getenv("CONFIG_ENERGY"));
 529     cout << "Energy set to " << energy << " GeV" << endl;
 530   }
 531
 532   // Random Number seed
 533   if (gSystem->Getenv("CONFIG_SEED"))
 534   {
 535     seed = atoi(gSystem->Getenv("CONFIG_SEED"));
 536   }
 537 }