Writing of FLUKA input cards for physics configuration is now the responsibility of
[u/mrichter/AliRoot.git] / TFluka / TFlukaConfigOption.cxx
CommitLineData
1df5fa54 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$*/
17
18#include "TFlukaConfigOption.h"
fb2cbbec 19#include "TFlukaMCGeometry.h"
20#include "TFluka.h"
21#include "TFlukaCerenkov.h"
22
23#include <TString.h>
24#include <TObjArray.h>
25#include <TVirtualMC.h>
26#include <TGeoMaterial.h>
27
28Float_t TFlukaConfigOption::fgMatMin(-1.);
29Float_t TFlukaConfigOption::fgMatMax(-1.);
30FILE* TFlukaConfigOption::fgFile(0x0);
31TFlukaMCGeometry* TFlukaConfigOption::fgGeom(0x0);
32
33Double_t TFlukaConfigOption::fgDCutValue[11];
34Int_t TFlukaConfigOption::fgDProcessFlag[15];
35
36
1df5fa54 37ClassImp(TFlukaConfigOption)
38
39
40TFlukaConfigOption::TFlukaConfigOption()
41{
42 // Default constructor
fb2cbbec 43 fMedium = -1;
44 fCMatMin = -1.;
45 fCMatMax = -1.;
46 Int_t i;
47 for (i = 0; i < 11; i++) fCutValue[i] = -1.;
48 for (i = 0; i < 15; i++) fProcessFlag[i] = -1;
1df5fa54 49}
50
51
fb2cbbec 52TFlukaConfigOption::TFlukaConfigOption(Int_t medium)
1df5fa54 53{
54 // Constructor
fb2cbbec 55 fMedium = medium;
56 fCMatMin = -1.;
57 fCMatMax = -1.;
58 Int_t i;
59 for (i = 0; i < 11; i++) fCutValue[i] = -1.;
60 for (i = 0; i < 15; i++) fProcessFlag[i] = -1;
1df5fa54 61}
62
fb2cbbec 63void TFlukaConfigOption::SetCut(const char* flagname, Double_t val)
1df5fa54 64{
fb2cbbec 65 // Set a cut value
66 const TString cuts[11] =
67 {"CUTGAM", "CUTELE", "CUTNEU", "CUTHAD", "CUTMUO", "BCUTE", "BCUTM", "DCUTE", "DCUTM", "PPCUTM", "TOFMAX"};
68 Int_t i;
69 for (i = 0; i < 11; i++) {
70 if (cuts[i].CompareTo(flagname) == 0) {
71 fCutValue[i] = val;
72 if (fMedium == -1) fgDCutValue[i] = val;
73 break;
74 }
75 }
1df5fa54 76}
77
fb2cbbec 78void TFlukaConfigOption::SetProcess(const char* flagname, Int_t flag)
79{
80 // Set a process flag
81 const TString process[15] =
82 {"DCAY", "PAIR", "COMP", "PHOT", "PFIS", "DRAY", "ANNI", "BREM", "MUNU", "CKOV",
83 "HADR", "LOSS", "MULS", "RAYL", "STRA"};
84 Int_t i;
85 for (i = 0; i < 15; i++) {
86 if (process[i].CompareTo(flagname) == 0) {
87 fProcessFlag[i] = flag;
88 if (fMedium == -1) fgDProcessFlag[i] = flag;
89 break;
90 }
91 }
92}
1df5fa54 93
fb2cbbec 94void TFlukaConfigOption::WriteFlukaInputCards()
1df5fa54 95{
fb2cbbec 96 // Write the FLUKA input cards for the set of process flags and cuts
97 //
98 //
99 if (fMedium > -1) {
100 fprintf(fgFile,"*\n*Material specific process and cut settings for #%8d \n", fMedium);
101 fCMatMin = fMedium;
102 fCMatMax = fMedium;
103 } else {
104 fprintf(fgFile,"*\n*Global process and cut settings \n");
105 fCMatMin = fgMatMin;
106 fCMatMax = fgMatMax;
107 }
108
109//
110// Handle Process Flags
111//
112 if (fProcessFlag[kDCAY] != -1) ProcessDCAY();
113 if (fProcessFlag[kPAIR] != -1) ProcessPAIR();
114 if (fProcessFlag[kBREM] != -1) ProcessBREM();
115 if (fProcessFlag[kCOMP] != -1) ProcessCOMP();
116 if (fProcessFlag[kPHOT] != -1) ProcessPHOT();
117 if (fProcessFlag[kPFIS] != -1) ProcessPFIS();
118 if (fProcessFlag[kANNI] != -1) ProcessANNI();
119 if (fProcessFlag[kMUNU] != -1) ProcessMUNU();
120 if (fProcessFlag[kHADR] != -1) ProcessHADR();
121 if (fProcessFlag[kMULS] != -1) ProcessMULS();
122 if (fProcessFlag[kRAYL] != -1) ProcessRAYL();
123
124 if (fProcessFlag[kLOSS] != -1 || fProcessFlag[kDRAY] != -1) ProcessLOSS();
125 if ((fMedium == -1 && fProcessFlag[kCKOV] > 0) || (fMedium > -1 && fProcessFlag[kCKOV] != -1)) ProcessCKOV();
126
127//
128// Handle Cuts
129//
130 if (fCutValue[kCUTGAM] >= 0.) ProcessCUTGAM();
131 if (fCutValue[kCUTELE] >= 0.) ProcessCUTELE();
132 if (fCutValue[kCUTNEU] >= 0.) ProcessCUTNEU();
133 if (fCutValue[kCUTHAD] >= 0.) ProcessCUTHAD();
134 if (fCutValue[kCUTMUO] >= 0.) ProcessCUTMUO();
135
136 if (fCutValue[kTOFMAX] >= 0.) ProcessTOFMAX();
1df5fa54 137}
138
fb2cbbec 139void TFlukaConfigOption::ProcessDCAY()
1df5fa54 140{
fb2cbbec 141 // Process DCAY option
142 fprintf(fgFile,"*\n* --- DCAY --- Decays. Flag = %5d\n", fProcessFlag[kDCAY]);
143 if (fProcessFlag[kDCAY] == 0) {
144 printf("Decays cannot be switched off \n");
145 } else {
146 fprintf(fgFile, "* Decays are on by default\n");
147 }
148}
149
150
151void TFlukaConfigOption::ProcessPAIR()
152{
153 // Process PAIR option
154 fprintf(fgFile,"*\n* --- PAIR --- Pair production by gammas, muons and hadrons. Flag = %5d, PPCUTM = %13.4g \n",
155 fProcessFlag[kPAIR], fCutValue[kPPCUTM]);
156 //
157 // gamma -> e+ e-
158 //
159 if (fProcessFlag[kPAIR] > 0) {
160 fprintf(fgFile,"EMFCUT %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fPHOT-THR\n",0., 0., 0., fCMatMin, fCMatMax, 1.);
161 } else {
162 fprintf(fgFile,"EMFCUT %10.1f%10.1f%10.4g%10.1f%10.1f%10.1fPHOT-THR\n",0., 0., 1e10, fCMatMin, fCMatMax, 1.);
163 }
164
165 //
166 // Direct pair production by Muons and Hadrons
167 //
168 //
169 // Attention ! This card interferes with BREM
170 //
171
172 if (fProcessFlag[kBREM] == -1 ) fProcessFlag[kBREM] = fgDProcessFlag[kBREM];
173 if (fCutValue[kBCUTM] == -1.) fCutValue[kBCUTM] = fgDCutValue[kBCUTM];
174
175
176 Float_t flag = -3.;
177 if (fProcessFlag[kPAIR] > 0 && fProcessFlag[kBREM] == 0) flag = 1.;
178 if (fProcessFlag[kPAIR] == 0 && fProcessFlag[kBREM] > 0) flag = 2.;
179 if (fProcessFlag[kPAIR] > 0 && fProcessFlag[kBREM] > 0) flag = 3.;
180 if (fProcessFlag[kPAIR] == 0 && fProcessFlag[kBREM] == 0) flag = -3.;
181 // Flag BREM card as handled
182 fProcessFlag[kBREM] = -1;
183
184 //
185 // Energy cut for pair prodution
186 //
187 Float_t cutP = fCutValue[kPPCUTM];
188 if (fCutValue[kPPCUTM] == -1.) cutP = fgDCutValue[kPPCUTM];
189 // In G3 this is the cut on the total energy of the e+e- pair
190 // In FLUKA the cut is on the kinetic energy of the electron and poistron
191 cutP = cutP / 2. - 0.51099906e-3;
192 if (cutP < 0.) cutP = 0.;
193 // No explicite generation of e+/e-
194 if (fProcessFlag[kPAIR] == 2) cutP = -1.;
195 //
196 // Energy cut for bremsstrahlung
197 //
198 Float_t cutB = 0.;
199 if (flag > 1.) {
200 fprintf(fgFile,"*\n* +++ BREM --- Bremsstrahlung by muons/hadrons. Flag = %5d, BCUTM = %13.4g \n",
201 fProcessFlag[kBREM], fCutValue[kBCUTM]);
202
203 cutB = fCutValue[kBCUTM];
204 // No explicite production of gammas
205 if (fProcessFlag[kBREM] == 2) cutB = -1.;
206 }
207
208 fprintf(fgFile,"PAIRBREM %10.1f%10.4g%10.4g%10.1f%10.1f\n",flag, cutP, cutB, fCMatMin, fCMatMax);
209}
210
211
212void TFlukaConfigOption::ProcessBREM()
213{
214 // Process BREM option
215 fprintf(fgFile,"*\n* --- BREM --- Bremsstrahlung by e+/- and muons/hadrons. Flag = %5d, BCUTE = %13.4g, BCUTM = %13.4g \n",
216 fProcessFlag[kBREM], fCutValue[kBCUTE], fCutValue[kBCUTM]);
217
218 //
219 // e+/- -> e+/- gamma
220 //
221 Float_t cutB = fCutValue[kBCUTE];
222 if (fCutValue[kBCUTE] == -1.) cutB = fgDCutValue[kBCUTE];
223
224
225 if (fProcessFlag[kBREM] > 0) {
226 fprintf(fgFile,"EMFCUT %10.4g%10.1f%10.1f%10.1f%10.1f%10.1fELPO-THR\n",cutB, 0., 0., fCMatMin, fCMatMax, 1.);
227 } else {
228 fprintf(fgFile,"EMFCUT %10.4g%10.1f%10.1f%10.1f%10.1f%10.1fELPO-THR\n",1.e10, 0., 0., fCMatMin, fCMatMax, 1.);
229 }
230
231 //
232 // Bremsstrahlung by muons and hadrons
233 //
234 cutB = fCutValue[kBCUTM];
235 if (fCutValue[kBCUTM] == -1.) cutB = fgDCutValue[kBCUTM];
236 if (fProcessFlag[kBREM] == 2) cutB = -1.;
237 Float_t flag = 2.;
238 if (fProcessFlag[kBREM] == 0) flag = -2.;
239
240 fprintf(fgFile,"PAIRBREM %10.1f%10.4g%10.4g%10.1f%10.1f\n", flag, 0., cutB, fCMatMin, fCMatMax);
241}
242
243void TFlukaConfigOption::ProcessCOMP()
244{
245 // Process COMP option
246 fprintf(fgFile,"*\n* --- COMP --- Compton scattering Flag = %5d \n", fProcessFlag[kCOMP]);
247
248 //
249 // Compton scattering
250 //
251
252 if (fProcessFlag[kCOMP] > 0) {
253 fprintf(fgFile,"EMFCUT %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fPHOT-THR\n",0. , 0., 0., fCMatMin, fCMatMax, 1.);
254 } else {
255 fprintf(fgFile,"EMFCUT %10.4g%10.1f%10.1f%10.1f%10.1f%10.1fPHOT-THR\n",1.e10, 0., 0., fCMatMin, fCMatMax, 1.);
256 }
257}
258
259void TFlukaConfigOption::ProcessPHOT()
260{
261 // Process PHOS option
262 fprintf(fgFile,"*\n* --- PHOT --- Photoelectric effect. Flag = %5d\n", fProcessFlag[kPHOT]);
263
264 //
265 // Photoelectric effect
266 //
267
268 if (fProcessFlag[kPHOT] > 0) {
269 fprintf(fgFile,"EMFCUT %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fPHOT-THR\n",0. , 0., 0., fCMatMin, fCMatMax, 1.);
270 } else {
271 fprintf(fgFile,"EMFCUT %10.1f%10.4g%10.1f%10.1f%10.1f%10.1fPHOT-THR\n",0., 1.e10, 0., fCMatMin, fCMatMax, 1.);
272 }
273}
274
275void TFlukaConfigOption::ProcessANNI()
276{
277 // Process ANNI option
278 fprintf(fgFile,"*\n* --- ANNI --- Positron annihilation. Flag = %5d \n", fProcessFlag[kANNI]);
279
280 //
281 // Positron annihilation
282 //
283
284 if (fProcessFlag[kANNI] > 0) {
285 fprintf(fgFile,"EMFCUT %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fANNH-THR\n",0. , 0., 0., fCMatMin, fCMatMax, 1.);
286 } else {
287 fprintf(fgFile,"EMFCUT %10.4g%10.1f%10.1f%10.1f%10.1f%10.1fANNH-THR\n",1.e10, 0., 0., fCMatMin, fCMatMax, 1.);
288 }
289}
290
291
292void TFlukaConfigOption::ProcessPFIS()
293{
294 // Process PFIS option
295 fprintf(fgFile,"*\n* --- PFIS --- Photonuclear interaction Flag = %5d\n", fProcessFlag[kPFIS]);
296
297 //
298 // Photonuclear interactions
299 //
300
301 if (fProcessFlag[kPFIS] > 0) {
302 fprintf(fgFile,"PHOTONUC %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n",(Float_t) fProcessFlag[kPFIS], 0., 0., fCMatMin, fCMatMax, 1.);
303 } else {
304 fprintf(fgFile,"PHOTONUC %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n",-1. , 0., 0., fCMatMin, fCMatMax, 1.);
305 }
306}
307
308void TFlukaConfigOption::ProcessMUNU()
309{
310 // Process MUNU option
311 fprintf(fgFile,"*\n* --- MUNU --- Muon nuclear interaction. Flag = %5d\n", fProcessFlag[kMUNU]);
312
313 //
314 // Muon nuclear interactions
315 //
316 if (fProcessFlag[kMUNU] > 0) {
317 fprintf(fgFile,"MUPHOTON %10.1f%10.3f%10.3f%10.1f%10.1f%10.1f\n",(Float_t )fProcessFlag[kMUNU], 0.25, 0.75, fCMatMin, fCMatMax, 1.);
318 } else {
319 fprintf(fgFile,"MUPHOTON %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n",-1. , 0., 0., fCMatMin, fCMatMax, 1.);
320 }
321}
322
323void TFlukaConfigOption::ProcessRAYL()
324{
325 // Process RAYL option
326 fprintf(fgFile,"*\n* --- RAYL --- Rayleigh Scattering. Flag = %5d\n", fProcessFlag[kRAYL]);
327
328 //
329 // Rayleigh scattering
330 //
331 Int_t nreg;
332 Int_t* reglist = fgGeom->GetMaterialList(fMedium, nreg);
333 //
334 // Loop over regions of a given material
335 for (Int_t k = 0; k < nreg; k++) {
336 Float_t ireg = reglist[k];
337 if (fProcessFlag[kRAYL] > 0) {
338 fprintf(fgFile,"EMFRAY %10.1f%10.1f%10.1f%10.1f\n", 1., ireg, ireg, 1.);
339 } else {
340 fprintf(fgFile,"EMFRAY %10.1f%10.1f%10.1f%10.1f\n", 3., ireg, ireg, 1.);
341 }
342 }
343}
344
345void TFlukaConfigOption::ProcessCKOV()
346{
347 // Process CKOV option
348 fprintf(fgFile,"*\n* --- CKOV --- Cerenkov Photon production. %5d\n", fProcessFlag[kCKOV]);
349
350 //
351 // Cerenkov photon production
352 //
353
354 TFluka* fluka = (TFluka*) gMC;
355 TObjArray *matList = fluka->GetFlukaMaterials();
356 Int_t nmaterial = matList->GetEntriesFast();
357 for (Int_t im = 0; im < nmaterial; im++)
358 {
359 TGeoMaterial* material = dynamic_cast<TGeoMaterial*> (matList->At(im));
360 Int_t idmat = material->GetIndex();
361//
362// Check if global option
363 if (fMedium != -1 && idmat != fMedium) continue;
364
365 TFlukaCerenkov* cerenkovProp;
366 if (!(cerenkovProp = dynamic_cast<TFlukaCerenkov*>(material->GetCerenkovProperties()))) continue;
367 //
368 // This medium has Cerenkov properties
369 //
370 //
371 if (fMedium == -1 || (fMedium != -1 && fProcessFlag[kCKOV] > 0)) {
372 // Write OPT-PROD card for each medium
373 Float_t emin = cerenkovProp->GetMinimumEnergy();
374 Float_t emax = cerenkovProp->GetMaximumEnergy();
375 fprintf(fgFile, "OPT-PROD %10.4g%10.4g%10.4g%10.4g%10.4g%10.4gCERENKOV\n", emin, emax, 0.,
376 Float_t(idmat), Float_t(idmat), 0.);
377 //
378 // Write OPT-PROP card for each medium
379 // Forcing FLUKA to call user routines (queffc.cxx, rflctv.cxx, rfrndx.cxx)
380 //
381 fprintf(fgFile, "OPT-PROP %10.4g%10.4g%10.4g%10.1f%10.1f%10.1fWV-LIMIT\n",
382 cerenkovProp->GetMinimumWavelength(), cerenkovProp->GetMaximumWavelength(), cerenkovProp->GetMaximumWavelength(),
383 Float_t(idmat), Float_t(idmat), 0.0);
384
385 if (cerenkovProp->IsMetal()) {
386 fprintf(fgFile, "OPT-PROP %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fMETAL\n", -100., -100., -100.,
387 Float_t(idmat), Float_t(idmat), 0.0);
388 } else {
389 fprintf(fgFile, "OPT-PROP %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f\n", -100., -100., -100.,
390 Float_t(idmat), Float_t(idmat), 0.0);
391 }
392
393
394 for (Int_t j = 0; j < 3; j++) {
395 fprintf(fgFile, "OPT-PROP %10.1f%10.1f%10.1f%10.1f%10.1f%10.1f&\n", -100., -100., -100.,
396 Float_t(idmat), Float_t(idmat), 0.0);
397 }
398 // Photon detection efficiency user defined
399
400 if (cerenkovProp->IsSensitive())
401 fprintf(fgFile, "OPT-PROP %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fSENSITIV\n", -100., -100., -100.,
402 Float_t(idmat), Float_t(idmat), 0.0);
403 } else {
404 fprintf(fgFile,"OPT-PROD %10.1f%10.1f%10.1f%10.1f%10.1f%10.1fCERE-OFF\n",0., 0., 0., fCMatMin, fCMatMax, 1.);
405 }
406 }
407}
408
409
410void TFlukaConfigOption::ProcessHADR()
411{
412 // Process HADR option
413 fprintf(fgFile,"*\n* --- HADR --- Hadronic interactions. Flag = %5d\n", fProcessFlag[kHADR]);
414
415 if (fProcessFlag[kHADR] > 0) {
416 fprintf(fgFile,"*\n*Hadronic interaction is ON by default in FLUKA\n");
417 } else {
418 if (fMedium != -1) printf("Hadronic interactions cannot be switched off material by material !\n");
419 fprintf(fgFile,"THRESHOL %10.1f%10.1f%10.1f%10.1e%10.1f\n",0., 0., 0., 1.e10, 0.);
420 }
421}
422
423
424
425void TFlukaConfigOption::ProcessMULS()
426{
427 // Process MULS option
428 fprintf(fgFile,"*\n* --- MULS --- Muliple Scattering. Flag = %5d\n", fProcessFlag[kMULS]);
429 //
430 // Multiple scattering
431 //
432 if (fProcessFlag[kMULS] > 0) {
433 fprintf(fgFile,"*\n*Multiple scattering is ON by default in FLUKA\n");
434 } else {
435 fprintf(fgFile,"MULSOPT %10.1f%10.1f%10.1f%10.1f%10.1f\n",0., 3., 3., fCMatMin, fCMatMax);
436 }
437}
438
439void TFlukaConfigOption::ProcessLOSS()
440{
441 // Process LOSS option
442 fprintf(fgFile,"*\n* --- LOSS --- Ionisation energy loss. Flags: LOSS = %5d, DRAY = %5d, STRA = %5d; Cuts: DCUTE = %13.4g, DCUTM = %13.4g \n",
443 fProcessFlag[kLOSS], fProcessFlag[kDRAY], fProcessFlag[kSTRA], fCutValue[kDCUTE], fCutValue[kDCUTM]);
444 //
445 // Ionisation energy loss
446 //
447 //
448 // Impose consistency
449
450 if (fProcessFlag[kLOSS] == 1 || fProcessFlag[kLOSS] == 3) {
451 fProcessFlag[kDRAY] = 1;
452 } else if (fProcessFlag[kLOSS] == 2) {
453 fProcessFlag[kDRAY] = 0;
454 fCutValue[kDCUTE] = 1.e10;
455 fCutValue[kDCUTM] = 1.e10;
456 } else {
457 if (fProcessFlag[kDRAY] == 1) {
458 fProcessFlag[kLOSS] = 1;
459 } else if (fProcessFlag[kDRAY] == 0) {
460 fProcessFlag[kLOSS] = 2;
461 fCutValue[kDCUTE] = 1.e10;
462 fCutValue[kDCUTM] = 1.e10;
463 }
464 }
465
466 if (fCutValue[kDCUTE] == -1.) fCutValue[kDCUTE] = fgDCutValue[kDCUTE];
467 if (fCutValue[kDCUTM] == -1.) fCutValue[kDCUTM] = fgDCutValue[kDCUTM];
468
469 Float_t cutM = fCutValue[kDCUTM];
470
471
472 if (fProcessFlag[kSTRA] == -1) fProcessFlag[kSTRA] = fgDProcessFlag[kSTRA];
473 if (fProcessFlag[kSTRA] == 0) fProcessFlag[kSTRA] = 1;
474 Float_t stra = (Float_t) fProcessFlag[kSTRA];
475
476
477 if (fProcessFlag[kLOSS] == 1 || fProcessFlag[kLOSS] == 3) {
478//
479// Restricted energy loss fluctuations
480//
481 fprintf(fgFile,"IONFLUCT %10.1f%10.1f%10.1f%10.1f%10.1f\n", 1., 1., stra, fCMatMin, fCMatMax);
482 fprintf(fgFile,"DELTARAY %10.4g%10.1f%10.1f%10.1f%10.1f%10.1f\n", cutM, 0., 0., fCMatMin, fCMatMax, 1.);
483 } else if (fProcessFlag[kLOSS] == 4) {
484//
485// No fluctuations
486//
487 fprintf(fgFile,"IONFLUCT %10.1f%10.1f%10.1f%10.1f%10.1f\n",-1., -1., stra, fCMatMin, fCMatMax);
488 fprintf(fgFile,"DELTARAY %10.4g%10.1f%10.1f%10.1f%10.1f%10.1f\n", 1.e10, 0., 0., fCMatMin, fCMatMax, 1.);
489 } else {
490//
491// Full fluctuations
492//
493 fprintf(fgFile,"IONFLUCT %10.1f%10.1f%10.1f%10.1f%10.1f\n",1., 1., stra, fCMatMin, fCMatMax);
494 fprintf(fgFile,"DELTARAY %10.4g%10.1f%10.1f%10.1f%10.1f%10.1f\n", 1.e10, 0., 0., fCMatMin, fCMatMax, 1.);
495 }
496}
497
498
499void TFlukaConfigOption::ProcessCUTGAM()
500{
501// Cut on gammas
502//
503 fprintf(fgFile,"*\n*Cut for Gammas. CUTGAM = %13.4g\n", fCutValue[kCUTGAM]);
504 if (fMedium == -1) {
505 fprintf(fgFile,"EMFCUT %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n",
506 0., fCutValue[kCUTGAM], 0., 0., Float_t(fgGeom->NofVolumes()), 1.);
507 } else {
508 Int_t nreg, *reglist;
509 Float_t ireg;
510 reglist = fgGeom->GetMaterialList(fMedium, nreg);
511 // Loop over regions of a given material
512 for (Int_t k = 0; k < nreg; k++) {
513 ireg = reglist[k];
514 fprintf(fgFile,"EMFCUT %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n", 0.,fCutValue[kCUTGAM], 0., ireg, ireg, 1.);
515 }
516 }
517}
518
519void TFlukaConfigOption::ProcessCUTELE()
520{
521// Cut on e+/e-
522//
523 fprintf(fgFile,"*\n*Cut for e+/e-. CUTELE = %13.4g\n", fCutValue[kCUTELE]);
524 if (fMedium == -1) {
525 fprintf(fgFile,"EMFCUT %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n",
526 -fCutValue[kCUTELE], 0., 0., 0., Float_t(fgGeom->NofVolumes()), 1.);
527 } else {
528 Int_t nreg, *reglist;
529 Float_t ireg;
530 reglist = fgGeom->GetMaterialList(fMedium, nreg);
531 // Loop over regions of a given material
532 for (Int_t k = 0; k < nreg; k++) {
533 ireg = reglist[k];
534 fprintf(fgFile,"EMFCUT %10.4g%10.4g%10.1f%10.1f%10.1f%10.1f\n", -fCutValue[kCUTELE], 0., 0., ireg, ireg, 1.);
535 }
536 }
537}
538
539void TFlukaConfigOption::ProcessCUTNEU()
540{
541 // Cut on neutral hadrons
542 fprintf(fgFile,"*\n*Cut for neutal hadrons. CUTNEU = %13.4g\n", fCutValue[kCUTNEU]);
543 if (fMedium == -1) {
544 Float_t cut = fCutValue[kCUTNEU];
545 // 8.0 = Neutron
546 // 9.0 = Antineutron
547 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 8.0, 9.0);
548 // 12.0 = Kaon zero long
549 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 12.0, 12.0);
550 // 17.0 = Lambda, 18.0 = Antilambda
551 // 19.0 = Kaon zero short
552 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 17.0, 19.0);
553 // 22.0 = Sigma zero, Pion zero, Kaon zero
554 // 25.0 = Antikaon zero
555 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 22.0, 25.0);
556 // 32.0 = Antisigma zero
557 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 32.0, 32.0);
558 // 34.0 = Xi zero
559 // 35.0 = AntiXi zero
560 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 34.0, 35.0);
561 // 47.0 = D zero
562 // 48.0 = AntiD zero
563 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 47.0, 48.0);
564 // 53.0 = Xi_c zero
565 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 53.0, 53.0);
566 // 55.0 = Xi'_c zero
567 // 56.0 = Omega_c zero
568 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 55.0, 56.0);
569 // 59.0 = AntiXi_c zero
570 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 59.0, 59.0);
571 // 61.0 = AntiXi'_c zero
572 // 62.0 = AntiOmega_c zero
573 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 61.0, 62.0);
574 } else {
575 printf("Cuts on neutral hadrons material by material not yet implemented !\n");
576 }
577}
578
579void TFlukaConfigOption::ProcessCUTHAD()
580{
581 // Cut on charged hadrons
582 fprintf(fgFile,"*\n*Cut for charge hadrons. CUTHAD = %13.4g\n", fCutValue[kCUTHAD]);
583 if (fMedium == -1) {
584 Float_t cut = fCutValue[kCUTHAD];
585 // 1.0 = Proton
586 // 2.0 = Antiproton
587 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 1.0, 2.0);
588 // 13.0 = Positive Pion, Negative Pion, Positive Kaon
589 // 16.0 = Negative Kaon
590 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 13.0, 16.0);
591 // 20.0 = Negative Sigma
592 // 21.0 = Positive Sigma
593 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 20.0, 21.0);
594 // 31.0 = Antisigma minus
595 // 33.0 = Antisigma plus
596 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 31.0, 31.0);
597 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 33.0, 33.0);
598 // 36.0 = Negative Xi, Positive Xi, Omega minus
599 // 39.0 = Antiomega
600 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 36.0, 39.0);
601 // 45.0 = D plus
602 // 46.0 = D minus
603 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 45.0, 46.0);
604 // 49.0 = D_s plus, D_s minus, Lambda_c plus
605 // 52.0 = Xi_c plus
606 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 49.0, 52.0);
607 // 54.0 = Xi'_c plus
608 // 60.0 = AntiXi'_c minus
609 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 54.0, 54.0);
610 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 60.0, 60.0);
611 // 57.0 = Antilambda_c minus
612 // 58.0 = AntiXi_c minus
613 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n", -cut, 57.0, 58.0);
614 } else {
615 printf("Cuts on charged hadrons material by material not yet implemented !\n");
616 }
617}
618
619void TFlukaConfigOption::ProcessCUTMUO()
620{
621 // Cut on muons
622 fprintf(fgFile,"*\n*Cut for muons. CUTMUO = %13.4g\n", fCutValue[kCUTMUO]);
623 Float_t cut = fCutValue[kCUTMUO];
624 if (fMedium == -1) {
625 fprintf(fgFile,"PART-THR %10.4g%10.1f%10.1f\n",-cut, 10.0, 11.0);
626 } else {
627 printf("Cuts on muons material by material not yet implemented !\n");
628 }
629
630
631}
632
633void TFlukaConfigOption::ProcessTOFMAX()
634{
635 // Cut time of flight
636 Float_t cut = fCutValue[kTOFMAX];
637 fprintf(fgFile,"*\n*Cut on time of flight. TOFMAX = %13.4g\n", fCutValue[kTOFMAX]);
638 fprintf(fgFile,"TIME-CUT %10.4g%10.1f%10.1f%10.1f%10.1f\n",cut*1.e9,0.,0.,-6.0,64.0);
1df5fa54 639}