All detectors active. Shish-kebab geometry for EMCAL
[u/mrichter/AliRoot.git] / 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   AliRunLoader* rl=0x0;
60                                                                                 
61   cout<<"Config.C: Creating Run Loader ..."<<endl;
62   rl = AliRunLoader::Open("galice.root",
63                             AliConfig::GetDefaultEventFolderName(),
64                             "recreate");
65   if (rl == 0x0)
66     {
67       gAlice->Fatal("Config.C","Can not instatiate the Run Loader");
68       return;
69     }
70   rl->SetCompressionLevel(2);
71   rl->SetNumberOfEventsPerFile(3);
72   gAlice->SetRunLoader(rl);
73                                                                                 
74   //
75   // Set External decayer
76   AliDecayer *decayer = new AliDecayerPythia();
77                                                                                
78   decayer->SetForceDecay(kAll);
79   decayer->Init();
80   gMC->SetExternalDecayer(decayer);
81
82   //
83   //
84   //
85   // Physics process control
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",1);
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   //=======================================================================
119   // ************* STEERING parameters FOR ALICE SIMULATION **************
120   // --- Specify event type to be tracked through the ALICE setup
121   // --- All positions are in cm, angles in degrees, and P and E in GeV
122   if (gSystem->Getenv("CONFIG_NPARTICLES"))
123       int     nParticles = atoi(gSystem->Getenv("CONFIG_NPARTICLES"));
124   else
125       int     nParticles = 1000;
126   
127   cout << "\t* Creating and configuring generator for " << nParticles 
128        << " particles..." << endl;
129   
130   AliGenHIJINGpara *gener = new AliGenHIJINGpara(nParticles);
131   
132   gener->SetMomentumRange(0., 999);
133   gener->SetPhiRange(0, 360);
134   // Set pseudorapidity range from -3 to 3.
135   Float_t thmin = EtaToTheta( 3.);   // theta min. <---> eta max
136   Float_t thmax = EtaToTheta(-3.);   // theta max. <---> eta min 
137   gener->SetThetaRange(thmin,thmax);
138   gener->SetOrigin(0, 0, 0);  //vertex position
139   gener->SetSigma(0, 0, 0);   //Sigma in (X,Y,Z) (cm) on IP position
140   gener->Init();
141   // 
142   // Activate this line if you want the vertex smearing to happen
143   // track by track
144   //
145
146
147      if (smag == k2kG) {
148         comment = comment.Append(" | L3 field 0.2 T");
149     } else if (smag == k4kG) {
150         comment = comment.Append(" | L3 field 0.4 T");
151     } else if (smag == k5kG) {
152         comment = comment.Append(" | L3 field 0.5 T");
153     }
154                                                                                 
155                                                                                 
156     if (srad == kGluonRadiation)
157     {
158         comment = comment.Append(" | Gluon Radiation On");
159                                                                                 
160     } else {
161         comment = comment.Append(" | Gluon Radiation Off");
162     }
163                                                                                 
164     if (sgeo == kHoles)
165     {
166         comment = comment.Append(" | Holes for PHOS/RICH");
167                                                                                 
168     } else {
169         comment = comment.Append(" | No holes for PHOS/RICH");
170     }
171                                                                                 
172     printf("\n \n Comment: %s \n \n", comment.Data());
173                                                                                 
174                                                                                 
175 // Field (L3 0.4 T)
176     AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., smag);
177     field->SetL3ConstField(0); //Using const. field in the barrel
178     rl->CdGAFile();
179     gAlice->SetField(field);
180  
181   Int_t   iABSO  = 1; 
182   Int_t   iCRT   = 1; 
183   Int_t   iDIPO  = 1; 
184   Int_t   iFMD   = 1; 
185   Int_t   iFRAME = 1; 
186   Int_t   iHALL  = 1; 
187   Int_t   iITS   = 1; 
188   Int_t   iMAG   = 1; 
189   Int_t   iMUON  = 1; 
190   Int_t   iPHOS  = 1; 
191   Int_t   iPIPE  = 1; 
192   Int_t   iPMD   = 1; 
193   Int_t   iRICH  = 1; 
194   Int_t   iSHIL  = 1; 
195   Int_t   iSTART = 1; 
196   Int_t   iTOF   = 1; 
197   Int_t   iTPC   = 1;
198   Int_t   iTRD   = 1; 
199   Int_t   iZDC   = 1; 
200   Int_t   iEMCAL = 0; 
201   Int_t   iVZERO = 1;
202  
203   cout << "\t* Creating the detectors ..." << endl;
204   //=================== Alice BODY parameters =============================
205     //=================== Alice BODY parameters =============================
206     AliBODY *BODY = new AliBODY("BODY", "Alice envelop");
207                                                                                 
208                                                                                 
209     if (iMAG)
210     {
211         //=================== MAG parameters ============================
212         // --- Start with Magnet since detector layouts may be depending ---
213         // --- on the selected Magnet dimensions ---
214         AliMAG *MAG = new AliMAG("MAG", "Magnet");
215     }
216                                                                                 
217                                                                                 
218     if (iABSO)
219     {
220         //=================== ABSO parameters ============================
221         AliABSO *ABSO = new AliABSOv0("ABSO", "Muon Absorber");
222     }
223                                                                                 
224     if (iDIPO)
225     {
226         //=================== DIPO parameters ============================
227                                                                                 
228         AliDIPO *DIPO = new AliDIPOv2("DIPO", "Dipole version 2");
229     }
230                                                                                 
231     if (iHALL)
232     {
233         //=================== HALL parameters ============================
234                                                                                 
235         AliHALL *HALL = new AliHALL("HALL", "Alice Hall");
236     }
237                                                                                 
238                                                                                 
239     if (iFRAME)
240     {
241         //=================== FRAME parameters ============================
242                                                                                 
243         AliFRAMEv2 *FRAME = new AliFRAMEv2("FRAME", "Space Frame");
244         if (sgeo == kHoles) {
245             FRAME->SetHoles(1);
246         } else {
247             FRAME->SetHoles(0);
248         }
249     }
250                                                                                 
251     if (iSHIL)
252     {
253         //=================== SHIL parameters ============================
254                                                                                 
255         AliSHIL *SHIL = new AliSHILv2("SHIL", "Shielding Version 2");
256     }
257                                                                                 
258                                                                                 
259     if (iPIPE)
260     {
261         //=================== PIPE parameters ============================
262                                                                                 
263         AliPIPE *PIPE = new AliPIPEv0("PIPE", "Beam Pipe");
264     }
265                                                                                 
266     if(iITS) {
267                                                                                 
268     //=================== ITS parameters ============================
269     //
270     // As the innermost detector in ALICE, the Inner Tracking System "impacts" on
271     // almost all other detectors. This involves the fact that the ITS geometry
272     // still has several options to be followed in parallel in order to determine
273     // the best set-up which minimizes the induced background. All the geometries
274     // available to date are described in the following. Read carefully the comments
275     // and use the default version (the only one uncommented) unless you are making
276     // comparisons and you know what you are doing. In this case just uncomment the
277     // ITS geometry you want to use and run Aliroot.
278     //
279     // Detailed geometries:
280     //
281     //
282     //AliITS *ITS  = new AliITSv5symm("ITS","Updated ITS TDR detailed version with symmetric services");
283     //
284     //AliITS *ITS  = new AliITSv5asymm("ITS","Updates ITS TDR detailed version with asymmetric services");
285     //
286         AliITSvPPRasymmFMD *ITS  = new AliITSvPPRasymmFMD("ITS","New ITS PPR detailed version with asymmetric services");
287         ITS->SetMinorVersion(2);  // don't touch this parameter if you're not an ITS developer
288         ITS->SetReadDet(kTRUE);   // don't touch this parameter if you're not an ITS developer
289     //    ITS->SetWriteDet("$ALICE_ROOT/ITS/ITSgeometry_vPPRasymm2.det");  // don't touch this parameter if you're not an ITS developer
290         ITS->SetThicknessDet1(200.);   // detector thickness on layer 1 must be in the range [100,300]
291         ITS->SetThicknessDet2(200.);   // detector thickness on layer 2 must be in the range [100,300]
292         ITS->SetThicknessChip1(200.);  // chip thickness on layer 1 must be in the range [150,300]
293         ITS->SetThicknessChip2(200.);  // chip thickness on layer 2 must be in the range [150,300]
294         ITS->SetRails(0);            // 1 --> rails in ; 0 --> rails out
295         ITS->SetCoolingFluid(1);   // 1 --> water ; 0 --> freon
296                                                                                 
297     // Coarse geometries (warning: no hits are produced with these coarse geometries and they unuseful
298     // for reconstruction !):
299     //
300     //
301     //AliITSvPPRcoarseasymm *ITS  = new AliITSvPPRcoarseasymm("ITS","New ITS PPR coarse version with asymmetric services");
302     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
303     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
304     //
305     //AliITS *ITS  = new AliITSvPPRcoarsesymm("ITS","New ITS PPR coarse version with symmetric services");
306     //ITS->SetRails(0);                // 1 --> rails in ; 0 --> rails out
307     //ITS->SetSupportMaterial(0);      // 0 --> Copper ; 1 --> Aluminum ; 2 --> Carbon
308     //
309     //
310     //
311     // Geant3 <-> EUCLID conversion
312     // ============================
313     //
314     // SetEUCLID is a flag to output (=1) or not to output (=0) both geometry and
315     // media to two ASCII files (called by default ITSgeometry.euc and
316     // ITSgeometry.tme) in a format understandable to the CAD system EUCLID.
317     // The default (=0) means that you dont want to use this facility.
318     //
319         ITS->SetEUCLID(0);
320     }
321                                                                                 
322     if (iTPC)
323     {
324         //============================ TPC parameters ================================
325         // --- This allows the user to specify sectors for the SLOW (TPC geometry 2)
326         // --- Simulator. SecAL (SecAU) <0 means that ALL lower (upper)
327         // --- sectors are specified, any value other than that requires at least one
328         // --- sector (lower or upper)to be specified!
329         // --- Reminder: sectors 1-24 are lower sectors (1-12 -> z>0, 13-24 -> z<0)
330         // ---           sectors 25-72 are the upper ones (25-48 -> z>0, 49-72 -> z<0)
331         // --- SecLows - number of lower sectors specified (up to 6)
332         // --- SecUps - number of upper sectors specified (up to 12)
333         // --- Sens - sensitive strips for the Slow Simulator !!!
334         // --- This does NOT work if all S or L-sectors are specified, i.e.
335         // --- if SecAL or SecAU < 0
336         //
337         //
338         //-----------------------------------------------------------------------------
339                                                                                 
340         //  gROOT->LoadMacro("SetTPCParam.C");
341         //  AliTPCParam *param = SetTPCParam();
342         AliTPC *TPC = new AliTPCv2("TPC", "Default");
343                                                                                 
344         // All sectors included
345         TPC->SetSecAL(-1);
346         TPC->SetSecAU(-1);
347                                                                                 
348     }
349                                                                                 
350                                                                                 
351     if (iTOF) {
352         //=================== TOF parameters ============================
353         AliTOF *TOF = new AliTOFv4T0("TOF", "normal TOF");
354     }
355                                                                                 
356                                                                                 
357     if (iRICH)
358     {
359         //=================== RICH parameters ===========================
360         AliRICH *RICH = new AliRICHv1("RICH", "normal RICH");
361                                                                                 
362     }
363                                                                                 
364                                                                                 
365     if (iZDC)
366     {
367         //=================== ZDC parameters ============================
368                                                                                 
369         AliZDC *ZDC = new AliZDCv2("ZDC", "normal ZDC");
370     }
371                                                                                 
372     if (iTRD)
373     {
374         //=================== TRD parameters ============================
375                                                                                 
376         AliTRD *TRD = new AliTRDv1("TRD", "TRD slow simulator");
377                                                                                 
378         // Select the gas mixture (0: 97% Xe + 3% isobutane, 1: 90% Xe + 10% CO2)
379         TRD->SetGasMix(1);
380         if (sgeo == kHoles) {
381             // With hole in front of PHOS
382             TRD->SetPHOShole();
383             // With hole in front of RICH
384             TRD->SetRICHhole();
385         }
386             // Switch on TR
387             AliTRDsim *TRDsim = TRD->CreateTR();
388     }
389                                                                                 
390     if (iFMD)
391     {
392         //=================== FMD parameters ============================
393         AliFMD *FMD = new AliFMDv1("FMD", "normal FMD");
394    }
395                                                                                 
396     if (iMUON)
397     {
398         //=================== MUON parameters ===========================
399                                                                                 
400         AliMUON *MUON = new AliMUONv1("MUON", "default");
401         MUON->AddGeometryBuilder(new AliMUONSt1GeometryBuilder(MUON));
402         MUON->AddGeometryBuilder(new AliMUONSt2GeometryBuilder(MUON));
403         MUON->AddGeometryBuilder(new AliMUONSlatGeometryBuilder(MUON));
404         MUON->AddGeometryBuilder(new AliMUONTriggerGeometryBuilder(MUON));
405     }
406     //=================== PHOS parameters ===========================
407                                                                                 
408     if (iPHOS)
409     {
410         AliPHOS *PHOS = new AliPHOSv1("PHOS", "IHEP");
411     }
412                                                                                 
413                                                                                 
414     if (iPMD)
415     {
416         //=================== PMD parameters ============================
417         AliPMD *PMD = new AliPMDv1("PMD", "normal PMD");
418     }
419                                                                                 
420     if (iSTART)
421     {
422         //=================== START parameters ============================
423         AliSTART *START = new AliSTARTv1("START", "START Detector");
424     }
425                                                                                 
426     if (iEMCAL)
427     {
428         //=================== EMCAL parameters ============================
429         AliEMCAL *EMCAL = new AliEMCALv2("EMCAL", "SHISH");
430     }
431                                                                                 
432      if (iCRT)
433     {
434         //=================== CRT parameters ============================
435         AliCRT *CRT = new AliCRTv0("CRT", "normal ACORDE");
436     }
437                                                                                 
438      if (iVZERO)
439     {
440         //=================== CRT parameters ============================
441         AliVZERO *VZERO = new AliVZEROv3("VZERO", "normal VZERO");
442     }
443                                                                                 
444                                                                                 
445 }
446                                                                                 
447 Float_t EtaToTheta(Float_t arg){
448   return (180./TMath::Pi())*2.*atan(exp(-arg));
449 }
450