TFluka/macro/FlukaConfig.C

   1 static Int_t    eventsPerRun = 100;
   2 enum PprGeo_t
   3 {
   4     kHoles, kNoHoles
   5 };
   6
   7 enum PprRad_t
   8 {
   9     kGluonRadiation, kNoGluonRadiation
  10 };
  11
  12 enum PprMag_t
  13 {
  14     k2kG, k4kG, k5kG
  15 };
  16
  17
  18 // This part for configuration
  19 static PprGeo_t sgeo = kHoles;
  20 static PprRad_t srad = kGluonRadiation;
  21 static PprMag_t smag = k5kG;
  22
  23 // Comment line
  24 static TString  comment;
  25
  26 // Functions
  27 Float_t EtaToTheta(Float_t arg);
  28
  29
  30 void Config()
  31 {
  32   cout << "==> Config.C..." << endl;
  33
  34   // Set Random Number seed
  35   gRandom->SetSeed(12345);
  36   cout<<"Seed for random number generation= "<<gRandom->GetSeed()<<endl;
  37
  38
  39
  40   // libraries required by fluka21
  41
  42   Bool_t isFluka = kTRUE;
  43   if (isFluka) {
  44     gSystem->Load("libGeom");
  45     cout << "\t* Loading TFluka..." << endl;
  46     gSystem->Load("libTFluka");
  47
  48     cout << "\t* Instantiating TFluka..." << endl;
  49     new  TFluka("C++ Interface to Fluka", 0/*verbositylevel*/);
  50   }
  51   else {
  52     cout << "\t* Loading Geant3..." << endl;
  53     gSystem->Load("libgeant321");
  54
  55     cout << "\t* Instantiating Geant3TGeo..." << endl;
  56     new     TGeant3TGeo("C++ Interface to Geant3");
  57   }
  58
  59   if(!AliCDBManager::Instance()->IsDefaultStorageSet()){
  60       AliCDBManager::Instance()->SetDefaultStorage("local://$ALICE_ROOT");
  61       AliCDBManager::Instance()->SetRun(0);
  62   }
  63
  64
  65   AliRunLoader* rl=0x0;
  66
  67   cout<<"Config.C: Creating Run Loader ..."<<endl;
  68   rl = AliRunLoader::Open("galice.root",
  69                             AliConfig::GetDefaultEventFolderName(),
  70                             "recreate");
  71   if (rl == 0x0)
  72     {
  73       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
  74       return;
  75     }
  76   rl->SetCompressionLevel(2);
  77   rl->SetNumberOfEventsPerFile(3);
  78   gAlice->SetRunLoader(rl);
  79
  80   //
  81   // Set External decayer
  82   AliDecayer *decayer = new AliDecayerPythia();
  83
  84   decayer->SetForceDecay(kAll);
  85   decayer->Init();
  86   gMC->SetExternalDecayer(decayer);
  87
  88   //
  89   //
  90   //
  91   // Physics process control
  92
  93   gMC->SetProcess("DCAY",1);
  94   gMC->SetProcess("PAIR",1);
  95   gMC->SetProcess("COMP",1);
  96   gMC->SetProcess("PHOT",1);
  97   gMC->SetProcess("PFIS",0);
  98   gMC->SetProcess("DRAY",1);
  99   gMC->SetProcess("ANNI",1);
 100   gMC->SetProcess("BREM",1);
 101   gMC->SetProcess("MUNU",1);
 102   gMC->SetProcess("CKOV",1);
 103   gMC->SetProcess("HADR",1);
 104   gMC->SetProcess("LOSS",2);
 105   gMC->SetProcess("MULS",1);
 106   gMC->SetProcess("RAYL",1);
 107
 108   Float_t cut = 1.e-3;        // 1MeV cut by default
 109   Float_t tofmax = 1.e10;
 110
 111   gMC->SetCut("CUTGAM", cut);
 112   gMC->SetCut("CUTELE", cut);
 113   gMC->SetCut("CUTNEU", cut);
 114   gMC->SetCut("CUTHAD", cut);
 115   gMC->SetCut("CUTMUO", cut);
 116   gMC->SetCut("BCUTE",  cut);
 117   gMC->SetCut("BCUTM",  cut);
 118   gMC->SetCut("DCUTE",  cut);
 119   gMC->SetCut("DCUTM",  cut);
 120   gMC->SetCut("PPCUTM", cut);
 121   gMC->SetCut("TOFMAX", tofmax);
 122
 123   //
 124   //=======================================================================
 125   // ************* STEERING parameters FOR ALICE SIMULATION **************
 126   // --- Specify event type to be tracked through the ALICE setup
 127   // --- All positions are in cm, angles in degrees, and P and E in GeV
 128   if (gSystem->Getenv("CONFIG_NPARTICLES"))
 129       int     nParticles = atoi(gSystem->Getenv("CONFIG_NPARTICLES"));
 130   else
 131       int     nParticles = 1000;
 132
 133   cout << "\t* Creating and configuring generator for " << nParticles
 134        << " particles..." << endl;
 135
 136   AliGenHIJINGpara *gener = new AliGenHIJINGpara(nParticles);
 137
 138   gener->SetMomentumRange(0., 999);
 139   gener->SetPhiRange(0, 360);
 140   // Set pseudorapidity range from -3 to 3.
 141   Float_t thmin = EtaToTheta( 3.);   // theta min. <---> eta max
 142   Float_t thmax = EtaToTheta(-3.);   // theta max. <---> eta min
 143   gener->SetThetaRange(thmin,thmax);
 144   gener->SetOrigin(0, 0, 0);  //vertex position
 145   gener->SetSigma(0, 0, 0);   //Sigma in (X,Y,Z) (cm) on IP position
 146   gener->Init();
 147   //
 148   // Activate this line if you want the vertex smearing to happen
 149   // track by track
 150   //
 151
 152
 153      if (smag == k2kG) {
 154         comment = comment.Append(" | L3 field 0.2 T");
 155     } else if (smag == k4kG) {
 156         comment = comment.Append(" | L3 field 0.4 T");
 157     } else if (smag == k5kG) {
 158         comment = comment.Append(" | L3 field 0.5 T");
 159     }
 160
 161
 162     if (srad == kGluonRadiation)
 163     {
 164         comment = comment.Append(" | Gluon Radiation On");
 165
 166     } else {
 167         comment = comment.Append(" | Gluon Radiation Off");
 168     }
 169
 170     if (sgeo == kHoles)
 171     {
 172         comment = comment.Append(" | Holes for PHOS/HMPID");
 173
 174     } else {
 175         comment = comment.Append(" | No holes for PHOS/HMPID");
 176     }
 177
 178     printf("\n \n Comment: %s \n \n", comment.Data());
 179
 180
 181 // Field (L3 0.4 T)
 182     AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., smag);
 183     field->SetL3ConstField(0); //Using const. field in the barrel
 184     rl->CdGAFile();
 185     gAlice->SetField(field);
 186
 187   Int_t   iABSO  = 1;
 188   Int_t   iACORDE   = 1;
 189   Int_t   iDIPO  = 1;
 190   Int_t   iFMD   = 1;
 191   Int_t   iFRAME = 1;
 192   Int_t   iHALL  = 1;
 193   Int_t   iITS   = 1;
 194   Int_t   iMAG   = 1;
 195   Int_t   iMUON  = 1;
 196   Int_t   iPHOS  = 1;
 197   Int_t   iPIPE  = 1;
 198   Int_t   iPMD   = 1;
 199   Int_t   iHMPID  = 1;
 200   Int_t   iSHIL  = 1;
 201   Int_t   iT0 = 1;
 202   Int_t   iTOF   = 1;
 203   Int_t   iTPC   = 1;
 204   Int_t   iTRD   = 1;
 205   Int_t   iZDC   = 1;
 206   Int_t   iEMCAL = 0;
 207   Int_t   iVZERO = 1;
 208
 209   cout << "\t* Creating the detectors ..." << endl;
 210   //=================== Alice BODY parameters =============================
 211     //=================== Alice BODY parameters =============================
 212     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
 213
 214
 215     if (iMAG)
 216     {
 217         //=================== MAG parameters ============================
 218         // --- Start with Magnet since detector layouts may be depending ---
 219         // --- on the selected Magnet dimensions ---
 220         AliMAG *MAG = new AliMAG("MAG", "Magnet");
 221     }
 222
 223
 224     if (iABSO)
 225     {
 226         //=================== ABSO parameters ============================
 227         AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
 228     }
 229
 230     if (iDIPO)
 231     {
 232         //=================== DIPO parameters ============================
 233
 234         AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
 235     }
 236
 237     if (iHALL)
 238     {
 239         //=================== HALL parameters ============================
 240
 241         AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
 242     }
 243
 244
 245     if (iFRAME)
 246     {
 247         //=================== FRAME parameters ============================
 248
 249         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
 250     }
 251
 252     if (iSHIL)
 253     {
 254         //=================== SHIL parameters ============================
 255
 256         AliSHIL *SHIL = new AliSHILv2("SHIL", "Shielding Version 2");
 257     }
 258
 259
 260     if (iPIPE)
 261     {
 262         //=================== PIPE parameters ============================
 263
 264         AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
 265     }
 266
 267     if(iITS) {
 268
 269     //=================== ITS parameters ============================
 270     //
 271     // As the innermost detector in ALICE, the Inner Tracking System "impacts" on
 272     // almost all other detectors. This involves the fact that the ITS geometry
 273     // still has several options to be followed in parallel in order to determine
 274     // the best set-up which minimizes the induced background. All the geometries
 275     // available to date are described in the following. Read carefully the comments
 276     // and use the default version (the only one uncommented) unless you are making
 277     // comparisons and you know what you are doing. In this case just uncomment the
 278     // ITS geometry you want to use and run Aliroot.
 279     //
 280     // Detailed geometries:
 281     //
 282     //
 283     //AliITS *ITS  = new AliITSv5symm("ITS","Updated ITS TDR detailed version with symmetric services");
 284     //
 285     //AliITS *ITS  = new AliITSv5asymm("ITS","Updates ITS TDR detailed version with asymmetric services");
 286     //
 287         AliITSvPPRasymmFMD *ITS  = new AliITSvPPRasymmFMD("ITS","New ITS PPR detailed version with asymmetric services");
 288         ITS->SetMinorVersion(2);  // don't touch this parameter if you're not an ITS developer
 289         ITS->SetReadDet(kFALSE);          // don't touch this parameter if you're not an ITS developer
 290     //    ITS->SetWriteDet("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymm2.det");  // don't touch this parameter if you're not an ITS developer
 291         ITS->SetThicknessDet1(200.);   // detector thickness on layer 1 must be in the range [100,300]
 292         ITS->SetThicknessDet2(200.);   // detector thickness on layer 2 must be in the range [100,300]
 293         ITS->SetThicknessChip1(150.);  // chip thickness on layer 1 must be in the range [150,300]
 294         ITS->SetThicknessChip2(150.);  // chip thickness on layer 2 must be in the range [150,300]
 295         ITS->SetRails(0);            // 1 --> rails in ; 0 --> rails out
 296         ITS->SetCoolingFluid(1);   // 1 --> water ; 0 --> freon
 297
 298     // Coarse geometries (warning: no hits are produced with these coarse geometries and they unuseful
 299     // for reconstruction !):
 300     //
 301     //
 302     //AliITSvPPRcoarseasymm *ITS  = new AliITSvPPRcoarseasymm("ITS","New ITS PPR coarse version with asymmetric services");
 303     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
 304     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
 305     //
 306     //AliITS *ITS  = new AliITSvPPRcoarsesymm("ITS","New ITS PPR coarse version with symmetric services");
 307     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
 308     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
 309     //
 310     //
 311     //
 312     // Geant3 <-> EUCLID conversion
 313     // ============================
 314     //
 315     // SetEUCLID is a flag to output (=1) or not to output (=0) both geometry and
 316     // media to two ASCII files (called by default ITSgeometry.euc and
 317     // ITSgeometry.tme) in a format understandable to the CAD system EUCLID.
 318     // The default (=0) means that you dont want to use this facility.
 319     //
 320         ITS->SetEUCLID(0);
 321     }
 322
 323     if (iTPC)
 324     {
 325       //============================ TPC parameters =====================
 326         AliTPC *TPC = new AliTPCv2("TPC", "Default");
 327     }
 328
 329
 330     if (iTOF) {
 331         //=================== TOF parameters ============================
 332         AliTOF *TOF = new AliTOFv5T0("TOF", "normal TOF");
 333         // Partial geometry: modules at 2,3,4,6,7,11,12,14,15,16
 334         // starting at 6h in positive direction
 335         //      Int_t TOFSectors[18]={-1,-1,0,0,0,-1,0,0,-1,-1,-1,0,0,-1,0,0,0,0};
 336         // Partial geometry: modules at 1,2,6,7,9,10,11,12,15,16,17
 337         // (ALICE numbering convention)
 338         Int_t TOFSectors[18]={-1,0,0,-1,-1,-1,0,0,-1,0,0,0,0,-1,-1,0,0,0};
 339         TOF->SetTOFSectors(TOFSectors);
 340     }
 341
 342
 343     if (iHMPID)
 344     {
 345         //=================== HMPID parameters ===========================
 346         AliHMPID *HMPID = new AliHMPIDv1("HMPID", "normal HMPID");
 347
 348     }
 349
 350
 351     if (iZDC)
 352     {
 353         //=================== ZDC parameters ============================
 354
 355         AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
 356     }
 357
 358     if (iTRD)
 359     {
 360         //=================== TRD parameters ============================
 361
 362         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
 363         AliTRDgeometry *geoTRD = TRD->GetGeometry();
 364         // Partial geometry: modules at 2,3,4,6,11,12,14,15
 365         // starting at 6h in positive direction
 366         geoTRD->SetSMstatus(0,0);
 367         geoTRD->SetSMstatus(1,0);
 368         geoTRD->SetSMstatus(5,0);
 369         geoTRD->SetSMstatus(7,0);
 370         geoTRD->SetSMstatus(8,0);
 371         geoTRD->SetSMstatus(9,0);
 372         geoTRD->SetSMstatus(10,0);
 373         geoTRD->SetSMstatus(13,0);
 374         geoTRD->SetSMstatus(16,0);
 375         geoTRD->SetSMstatus(17,0);
 376     }
 377
 378     if (iFMD)
 379     {
 380         //=================== FMD parameters ============================
 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         AliMUON *MUON = new AliMUONv1("MUON", "default");
 389     }
 390     //=================== PHOS parameters ===========================
 391
 392     if (iPHOS)
 393     {
 394         AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
 395     }
 396
 397
 398     if (iPMD)
 399     {
 400         //=================== PMD parameters ============================
 401         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
 402     }
 403
 404     if (iT0)
 405     {
 406         //=================== T0 parameters ============================
 407         AliT0 *T0 = new AliT0v1("T0", "T0 Detector");
 408     }
 409
 410     if (iEMCAL)
 411     {
 412         //=================== EMCAL parameters ============================
 413         AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "SHISH_77_TRD1_2X2_FINAL_110DEG");
 414     }
 415
 416      if (iACORDE)
 417     {
 418         //=================== ACORDE parameters ============================
 419         AliACORDE *ACORDE = new AliACORDEv0("ACORDE", "normal ACORDE");
 420     }
 421
 422      if (iVZERO)
 423     {
 424         //=================== ACORDE parameters ============================
 425         AliVZERO *VZERO = new AliVZEROv7("VZERO", "normal VZERO");
 426     }
 427 }
 428
 429 Float_t EtaToTheta(Float_t arg){
 430   return (180./TMath::Pi())*2.*atan(exp(-arg));
 431 }
 432