macros/Config_gener.C

   1 #if !defined(__CINT__) || defined(__MAKECINT__)
   2 #include <Riostream.h>
   3 #include <TRandom.h>
   4 #include <TSystem.h>
   5 #include <TVirtualMC.h>
   6 #include <TGeant3.h>
   7 #include "STEER/AliRunLoader.h"
   8 #include "STEER/AliRun.h"
   9 #include "STEER/AliConfig.h"
  10 #include "PYTHIA6/AliDecayerPythia.h"
  11 #include "EVGEN/AliGenCocktail.h"
  12 #include "EVGEN/AliGenFixed.h"
  13 #include "EVGEN/AliGenBox.h"
  14 #include "EVGEN/AliGenScan.h"
  15 #include "EVGEN/AliGenHIJINGpara.h"
  16 #include "THijing/AliGenHijing.h"
  17 #include "PYTHIA6/AliGenPythia.h"
  18 #include "THerwig/AliGenHerwig.h"
  19 #include "EVGEN/AliGenParam.h"
  20 #include "EVGEN/AliGenMUONlib.h"
  21 #include "EVGEN/AliGenPHOSlib.h"
  22 #include "EVGEN/AliGenGSIlib.h"
  23 #include "EVGEN/AliGenFLUKAsource.h"
  24 #include "EVGEN/AliGenExtFile.h"
  25 #include "EVGEN/AliGenHalo.h"
  26 #include "EVGEN/AliGenReaderTreeK.h"
  27 #include "STRUCT/AliBODY.h"
  28 #include "STRUCT/AliMAG.h"
  29 #endif
  30
  31 enum gentype_t {hijing, hijingParam, gun, box, pythia, herwig,
  32                 param1, param2, param3, param4,
  33                 cocktail, fluka, halo, ntuple, scan};
  34
  35 gentype_t gentype = herwig;
  36
  37 Int_t ntracks=1;
  38
  39 void Config()
  40 {
  41
  42   // Set Random Number seed
  43   gRandom->SetSeed(12345); //Set 0 to use the current time
  44   cout<<"Seed for random number generation= "<<gRandom->GetSeed()<<endl;
  45
  46
  47   // libraries required by geant321
  48 #if defined(__CINT__)
  49   gSystem->Load("libgeant321");
  50 #endif
  51
  52   new TGeant3TGeo("C++ Interface to Geant3");
  53
  54 //=======================================================================
  55 //  Create the output file
  56
  57   AliRunLoader* rl=0x0;
  58
  59   cout<<"Config.C: Creating Run Loader ..."<<endl;
  60   rl = AliRunLoader::Open("galice.root",
  61                           AliConfig::GetDefaultEventFolderName(),
  62                           "recreate");
  63   if (rl == 0x0)
  64     {
  65       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
  66       return;
  67     }
  68   rl->SetCompressionLevel(2);
  69   rl->SetNumberOfEventsPerFile(3);
  70   gAlice->SetRunLoader(rl);
  71
  72 //
  73 // Set External decayer
  74   AliDecayer* decayer = new AliDecayerPythia();
  75   decayer->SetForceDecay(kAll);
  76   decayer->Init();
  77   gMC->SetExternalDecayer(decayer);
  78
  79
  80 //
  81 //=======================================================================
  82 // ************* STEERING parameters FOR ALICE SIMULATION **************
  83 // --- Specify event type to be tracked through the ALICE setup
  84 // --- All positions are in cm, angles in degrees, and P and E in GeV
  85
  86
  87   gMC->SetProcess("DCAY",1);
  88   gMC->SetProcess("PAIR",1);
  89   gMC->SetProcess("COMP",1);
  90   gMC->SetProcess("PHOT",1);
  91   gMC->SetProcess("PFIS",0);
  92   gMC->SetProcess("DRAY",0);
  93   gMC->SetProcess("ANNI",1);
  94   gMC->SetProcess("BREM",1);
  95   gMC->SetProcess("MUNU",1);
  96   gMC->SetProcess("CKOV",1);
  97   gMC->SetProcess("HADR",1);
  98   gMC->SetProcess("LOSS",2);
  99   gMC->SetProcess("MULS",1);
 100   gMC->SetProcess("RAYL",1);
 101
 102   Float_t cut = 1.e-3;        // 1MeV cut by default
 103   Float_t tofmax = 1.e10;
 104
 105   gMC->SetCut("CUTGAM", cut);
 106   gMC->SetCut("CUTELE", cut);
 107   gMC->SetCut("CUTNEU", cut);
 108   gMC->SetCut("CUTHAD", cut);
 109   gMC->SetCut("CUTMUO", cut);
 110   gMC->SetCut("BCUTE",  cut);
 111   gMC->SetCut("BCUTM",  cut);
 112   gMC->SetCut("DCUTE",  cut);
 113   gMC->SetCut("DCUTM",  cut);
 114   gMC->SetCut("PPCUTM", cut);
 115   gMC->SetCut("TOFMAX", tofmax);
 116
 117
 118   AliGenerator * gGener = 0x0;
 119   switch(gentype)
 120     {
 121     case gun:
 122 //*********************************************
 123 // Example for Fixed Particle Gun
 124 //*********************************************
 125       {
 126         AliGenFixed *gener = new AliGenFixed(ntracks);
 127         gener->SetMomentum(50);
 128         gener->SetPhi(180.);
 129         gener->SetTheta(5.);
 130         gener->SetOrigin(0,0,0);        //vertex position
 131         gener->SetPart(13);             //GEANT particle type
 132         gGener = gener;
 133       }
 134       break;
 135     case box:
 136 //*********************************************
 137 // Example for Moving Particle Gun            *
 138 //*********************************************
 139       {
 140         AliGenBox *gener = new AliGenBox(ntracks);
 141         gener->SetMomentumRange(3,4);
 142         gener->SetPhiRange(0,360);
 143         gener->SetThetaRange(90, 180. );
 144         gener->SetOrigin(0,0,0);
 145         //vertex position
 146         gener->SetSigma(0,0,0);         //Sigma in (X,Y,Z) (cm) on IP position
 147         gener->SetPart(5);              //GEANT particle type
 148         gGener = gener;
 149       }
 150       break;
 151     case scan:
 152 //*********************************************
 153 // Scanning on a grid                         *
 154 //*********************************************
 155       {
 156         AliGenScan *gener = new AliGenScan(-1);
 157         gener->SetMomentumRange(4,4);
 158         gener->SetPhiRange(0,360);
 159         gener->SetThetaRange(9,9);
 160         //vertex position
 161         gener->SetSigma(6,6,0);         //Sigma in (X,Y,Z) (cm) on IP position
 162         gener->SetPart(5);
 163         gener->SetRange(20, -100, 100, 20, -100, 100, 1, 500, 500);
 164         gGener = gener;
 165       }
 166       break;
 167
 168     case hijingParam:
 169       {
 170         AliGenHIJINGpara *gener = new AliGenHIJINGpara(ntracks);
 171         gener->SetMomentumRange(0,999);
 172         gener->SetPhiRange(0,360);
 173         gener->SetThetaRange(2,10);
 174         gener->SetOrigin(0,0,0);        //vertex position
 175         gener->SetSigma(0,0,0);         //Sigma in (X,Y,Z) (cm) on IP position
 176         gGener = gener;
 177       }
 178       break;
 179     case hijing:
 180       {
 181         AliGenHijing *gener = new AliGenHijing(-1);
 182 // centre of mass energy
 183         gener->SetEnergyCMS(5500);
 184 // reference frame
 185         gener->SetReferenceFrame("CMS     ");
 186 // projectile
 187         gener->SetProjectile("A       ", 208, 82);
 188         gener->SetTarget    ("A       ", 208, 82);
 189 // impact parameter range
 190         gener->SetImpactParameterRange(0, 3.);
 191 // evaluate cross section before run
 192         gener->SetEvaluate(0);
 193 // tell hijing to keep the full parent child chain
 194         gener->KeepFullEvent();
 195 // enable jet quenching
 196         gener->SetJetQuenching(1);
 197 // enable shadowing
 198         gener->SetShadowing(1);
 199 // neutral pion and heavy particle decays switched off
 200         gener->SetDecaysOff(1);
 201 // trigger
 202         gener->SetTrigger(0);
 203 // kinematic selection
 204         gener->SetSelectAll(0);
 205 // momentum range
 206         gener->SetMomentumRange(0,999);
 207 // phi range
 208         gener->SetPhiRange(0.,360.);
 209 // theta range
 210         gener->SetThetaRange(0,180.);
 211 // select flavor (0: no, 4: charm+beauty, 5:beauty)
 212         gener->SetFlavor(0);
 213 //
 214         gener->SetOrigin(0., 0.0 ,0);
 215         gener->SetSigma(0,0,5.3);
 216         gener->SetVertexSmear(kPerEvent);
 217 // no tracking
 218         gener->SetTrackingFlag(0);
 219         gGener = gener;
 220       }
 221       break;
 222
 223     case pythia:
 224 //********************************************
 225 // Example for Charm  Production with Pythia *
 226 //********************************************
 227       {
 228         AliGenPythia *gener = new AliGenPythia(-1);
 229 //   final state kinematic cuts
 230         gener->SetMomentumRange(0,999);
 231         gener->SetPhiRange(0. ,360.);
 232         gener->SetThetaRange(0., 180.);
 233         gener->SetYRange(-10,10);
 234         gener->SetPtRange(0,100);
 235 //   vertex position and smearing
 236         gener->SetOrigin(0,0,0);       // vertex position
 237         gener->SetVertexSmear(kPerEvent);
 238         gener->SetSigma(0,0,5.6);      // Sigma in (X,Y,Z) (cm) on IP position
 239 //   Structure function. See the list in EVGEN/AliStructFuncType.h
 240         gener->SetStrucFunc(kGRVHO);
 241 // Select corection for nuclear structure functions
 242 //     gener->SetNuclei(208,208);
 243 //
 244 //   Process type. See the list in PYTHIA6/AliPythia.h
 245         gener->SetProcess(kPyBeauty);
 246 //
 247 //   Pt transfer of the hard scattering
 248         gener->SetPtHard(0.,5.);
 249 //   Decay type (semielectronic, semimuonic, nodecay)
 250         gener->SetForceDecay(kSemiElectronic);
 251 //   Centre of mass energy
 252         gener->SetEnergyCMS(5500.);
 253 //   No Tracking
 254         gener->SetTrackingFlag(0);
 255         gGener = gener;
 256       }
 257       break;
 258
 259     case herwig:
 260 //********************************************
 261 // Example for Charm  Production with Pythia *
 262 //********************************************
 263       {
 264         AliGenHerwig *gener = new AliGenHerwig(-1);
 265 //   final state kinematic cuts
 266         gener->SetMomentumRange(0,7000);
 267         gener->SetPhiRange(0. ,360.);
 268         gener->SetThetaRange(0., 180.);
 269         gener->SetYRange(-10,10);
 270         gener->SetPtRange(0,7000);
 271 //   vertex position and smearing
 272         gener->SetOrigin(0,0,0);       // vertex position
 273         gener->SetVertexSmear(kPerEvent);
 274         gener->SetSigma(0,0,5.6);      // Sigma in (X,Y,Z) (cm) on IP position
 275 //   Beam momenta
 276         gener->SetBeamMomenta(7000,7000);
 277 //   Beams
 278         gener->SetProjectile("P");
 279         gener->SetTarget("P");
 280 //   Structure function
 281         gener->SetStrucFunc(kGRVHO);
 282 //   Hard scatering
 283         gener->SetPtHardMin(200);
 284         gener->SetPtRMS(20);
 285 //   Min bias
 286         gener->SetProcess(8000);
 287 //   No Tracking
 288         gener->SetTrackingFlag(0);
 289         gGener = gener;
 290       }
 291       break;
 292
 293     case param1:
 294 //*******************************************************
 295 // Example for J/psi  Production from  Parameterisation
 296 // using default library (AliMUONlib)
 297 //*******************************************************
 298       {
 299         AliGenParam *gener =
 300           new AliGenParam(ntracks, AliGenMUONlib::kUpsilon);
 301         gener->SetMomentumRange(0,999);
 302         gener->SetPtRange(0,999);
 303         gener->SetPhiRange(0. , 360.);
 304         gener->SetYRange(2.5,4);
 305         gener->SetCutOnChild(1);
 306         gener->SetChildThetaRange(2,9);
 307         gener->SetOrigin(0,0,0);        //vertex position
 308         gener->SetSigma(0,0,5.3);       //Sigma in (X,Y,Z) (cm) on IP position
 309         gener->SetForceDecay(kDiMuon);
 310         gener->SetTrackingFlag(0);
 311         gGener = gener;
 312       }
 313       break;
 314
 315     case param2:
 316 //*******************************************************
 317 // Example for Omega  Production from  Parameterisation
 318 // specifying library.
 319 //*******************************************************
 320       {
 321         AliGenParam *gener = new AliGenParam(1000,new AliGenPHOSlib(),
 322                                              AliGenPHOSlib::kOmega);
 323         gener->SetWeighting(kNonAnalog);
 324         gener->SetForceDecay(kNoDecay);
 325         gener->SetPtRange(0,100);
 326         gener->SetThetaRange(45,135);
 327         gener->SetTrackingFlag(0);
 328         gGener = gener;
 329       }
 330       break;
 331
 332     case param3:
 333 //*******************************************************
 334 // Example for Upsilon  Production from  Parameterisation
 335 // specifying library.
 336 // GSI style
 337 //*******************************************************
 338       {
 339         AliGenParam *gener = new AliGenParam(1000,new AliGenGSIlib(),
 340                                              AliGenGSIlib::kUpsilon, "MUON");
 341         gener->SetMomentumRange(0,999);
 342         gener->SetPtRange(0,999);
 343         gener->SetPhiRange(0., 360.);
 344         gener->SetYRange(2.5,4);
 345         gener->SetCutOnChild(1);
 346         gener->SetChildThetaRange(2,9);
 347         gener->SetOrigin(0,0,0);        //vertex position
 348         gener->SetSigma(0,0,5.3);       //Sigma in (X,Y,Z) (cm) on IP position
 349         gener->SetForceDecay(kDiMuon);
 350         gener->SetTrackingFlag(0);
 351         gGener = gener;
 352       }
 353       break;
 354
 355     case param4:
 356 //*******************************************************
 357 // Example for Omega  Production from  Parameterisation
 358 // specifying library.
 359 // The alternative way.
 360 //*******************************************************
 361       {
 362         AliGenLib* Lib=new AliGenPHOSlib();
 363         Int_t iOmega = AliGenPHOSlib::kOmega;
 364         AliGenParam *gener = new AliGenParam(50, iOmega,
 365                                              Lib->GetPt(iOmega, ""),
 366                                              Lib->GetY (iOmega, ""),
 367                                              Lib->GetIp(iOmega, ""));
 368         gener->SetPtRange(0,999);
 369         gener->SetWeighting(kNonAnalog);
 370         gener->SetForceDecay(kNoDecay);
 371         gener->SetTrackingFlag(0);
 372         gGener = gener;
 373       }
 374       break;
 375
 376     case fluka:
 377 //*******************************************************
 378 // Example for a FLUKA Boundary Source                  *
 379 //*******************************************************
 380       {
 381         AliGenFLUKAsource *gener = new AliGenFLUKAsource(-1);
 382         gener->SetFileName("$(ALICE_ROOT)/data/all32.root");
 383         gener->SetPartFlag(9);
 384         gener->SetAgeMax(1.e-5);
 385 //  31.7 events
 386         gener->SetFraction(0.0315);
 387 //     gener->SetFraction(0.75*0.0315);
 388         rl->CdGAFile();
 389 //     gener->SetPartFlag(10);
 390         gener->SetMomentumRange(0,999);
 391         gener->SetPhiRange(0.,360.);
 392         gener->SetThetaRange(0., 180.);
 393         gener->SetAgeMax(1.e-5);
 394
 395 //  31.7 events
 396 //     gener->SetFraction(0.0315);
 397         gGener = gener;
 398       }
 399       break;
 400
 401     case ntuple:
 402 //*******************************************************
 403 // Example for reading from a external file                  *
 404 //*******************************************************
 405       {
 406         AliGenExtFile *gener = new AliGenExtFile(-1);
 407         gener->SetVertexSmear(kPerEvent);
 408         gener->SetTrackingFlag(1);
 409
 410         AliGenReaderTreeK * reader = new AliGenReaderTreeK();
 411         reader->SetFileName("$(ALICE_ROOT)/data/dtujet93.root");
 412         gener->SetReader(reader);
 413         gGener = gener;
 414       }
 415       break;
 416
 417     case halo:
 418 //*******************************************************
 419 // Example for Tunnel Halo Source                       *
 420 //*******************************************************
 421       {
 422         AliGenHalo *gener = new AliGenHalo(ntracks);
 423         gener->SetFileName("/h1/morsch/marsip/marsip5.mu");
 424         gGener = gener;
 425       }
 426       break;
 427
 428     case cocktail:
 429 //*******************************************************
 430 // Example for a Cocktail                               *
 431 //*******************************************************
 432       {
 433         AliGenCocktail *gener = new AliGenCocktail();
 434
 435         gener->SetPhiRange(0,360);
 436         gener->SetYRange(2.5,4);
 437         gener->SetThetaRange(2,9);
 438         gener->SetPtRange(0,10);
 439         gener->SetOrigin(0,0,0);        //vertex position
 440         gener->SetSigma(0,0,0);         //Sigma in (X,Y,Z) (cm) on IP position
 441         gener->SetMomentumRange(0,999);
 442
 443         AliGenParam *jpsi = new AliGenParam(1,AliGenMUONlib::kJpsi);
 444         jpsi->SetForceDecay(kDiMuon);
 445         jpsi->SetCutOnChild(1);
 446
 447
 448         AliGenFLUKAsource *bg = new AliGenFLUKAsource(-1);
 449         bg->AddFile("$(ALICE_ROOT)/data/all32.root");
 450         rl->CdGAFile();
 451         bg->SetPartFlag(9);
 452         bg->SetAgeMax(1.e-5);
 453 //  31.7 events
 454 //     gener->SetFraction(0.0315);
 455         bg->SetFraction(0.01*0.0315);
 456
 457         gener->AddGenerator(jpsi,"J/Psi", 1);
 458         gener->AddGenerator(bg,"Background",1);
 459
 460         gGener = gener;
 461       }
 462       break;
 463     }
 464
 465 // Activate this line if you want the vertex smearing to happen
 466 // track by track
 467 //
 468 // gener->SetVertexSmear(kPerTrack);
 469
 470   gGener->Init();
 471
 472   gAlice->SetField(-999,2);    //Specify maximum magnetic field in Tesla (neg. ==> default field)
 473
 474   Int_t iMAG=1;
 475   rl->CdGAFile();
 476
 477 //=================== Alice BODY parameters =============================
 478   AliBODY *BODY = new AliBODY("BODY","Alice envelop");
 479
 480
 481   if(iMAG) {
 482 //=================== MAG parameters ============================
 483 // --- Start with Magnet since detector layouts may be depending ---
 484 // --- on the selected Magnet dimensions ---
 485     AliMAG *MAG  = new AliMAG("MAG","Magnet");
 486   }
 487 }