03ca248b |
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 | /* |
17 | $Log$ |
b0d8df96 |
18 | Revision 1.12 2003/01/31 12:18:53 morsch |
19 | Corrected indices. (E. Futo) |
20 | |
c230803a |
21 | Revision 1.9 2002/12/06 12:41:29 morsch |
22 | Mess from last merge cleaned up. |
23 | |
a7e55c80 |
24 | Revision 1.8 2002/12/06 12:28:44 morsch |
25 | Region to media mapping corrected and improved. |
26 | |
6d4d27f2 |
27 | Revision 1.7 2002/12/06 12:21:32 morsch |
28 | User stepping methods added (E. Futo) |
29 | |
fa3d1cc7 |
30 | Revision 1.6 2002/11/21 18:40:06 iglez2 |
31 | Media handling added |
32 | |
27b2f7fe |
33 | Revision 1.5 2002/11/07 17:59:10 iglez2 |
34 | Included the geometry through geant4_vmc/FLUGG |
35 | |
bf3aa28e |
36 | Revision 1.4 2002/11/04 16:00:46 iglez2 |
37 | The conversion between ID and PDG now uses Fluka routines and arrays which is more consistent. |
38 | |
f9cb2fec |
39 | Revision 1.3 2002/10/22 15:12:14 alibrary |
40 | Introducing Riostream.h |
41 | |
eae0fe66 |
42 | Revision 1.2 2002/10/14 14:57:40 hristov |
43 | Merging the VirtualMC branch to the main development branch (HEAD) |
44 | |
b9d0a01d |
45 | Revision 1.1.2.8 2002/10/08 16:33:17 iglez2 |
46 | LSOUIT is set to true before the second call to flukam. |
47 | |
48 | Revision 1.1.2.7 2002/10/08 09:30:37 iglez2 |
49 | Solved stupid missing ; |
50 | |
51 | Revision 1.1.2.6 2002/10/07 13:40:22 iglez2 |
52 | First implementations of the PDG <--> Fluka Id conversion routines |
53 | |
54 | Revision 1.1.2.5 2002/09/26 16:26:03 iglez2 |
55 | Added verbosity |
56 | Call to gAlice->Generator()->Generate() |
57 | |
58 | Revision 1.1.2.4 2002/09/26 13:22:23 iglez2 |
59 | Naive implementation of ProcessRun and ProcessEvent |
60 | Opening/Closing of input file (fInputFileName) with FORTRAN unit 5 before/after the first call to flukam inside Init() |
61 | |
62 | Revision 1.1.2.3 2002/09/20 15:35:51 iglez2 |
63 | Modification of LFDRTR. Value is passed to FLUKA !!! |
64 | |
65 | Revision 1.1.2.2 2002/09/18 14:34:44 iglez2 |
66 | Revised version with all pure virtual methods implemented |
67 | |
68 | Revision 1.1.2.1 2002/07/24 08:49:41 alibrary |
69 | Adding TFluka to VirtualMC |
70 | |
71 | Revision 1.1 2002/07/05 13:10:07 morsch |
72 | First commit of Fluka interface. |
73 | |
03ca248b |
74 | */ |
75 | |
eae0fe66 |
76 | #include <Riostream.h> |
b9d0a01d |
77 | |
6d4d27f2 |
78 | #include "TClonesArray.h" |
03ca248b |
79 | #include "TFluka.h" |
b9d0a01d |
80 | #include "TCallf77.h" //For the fortran calls |
81 | #include "Fdblprc.h" //(DBLPRC) fluka common |
b9d0a01d |
82 | #include "Fepisor.h" //(EPISOR) fluka common |
fa3d1cc7 |
83 | #include "Ffinuc.h" //(FINUC) fluka common |
84 | #include "Fiounit.h" //(IOUNIT) fluka common |
85 | #include "Fpaprop.h" //(PAPROP) fluka common |
f9cb2fec |
86 | #include "Fpart.h" //(PART) fluka common |
fa3d1cc7 |
87 | #include "Ftrackr.h" //(TRACKR) fluka common |
6d4d27f2 |
88 | #include "Fpaprop.h" //(PAPROP) fluka common |
fa3d1cc7 |
89 | #include "Ffheavy.h" //(FHEAVY) fluka common |
b9d0a01d |
90 | |
fa3d1cc7 |
91 | #include "TVirtualMC.h" |
bf3aa28e |
92 | #include "TG4GeometryManager.h" //For the geometry management |
93 | #include "TG4DetConstruction.h" //For the detector construction |
94 | |
95 | #include "FGeometryInit.hh" |
fa3d1cc7 |
96 | #include "TLorentzVector.h" |
6d4d27f2 |
97 | #include "FlukaVolume.h" |
bf3aa28e |
98 | |
b9d0a01d |
99 | // Fluka methods that may be needed. |
100 | #ifndef WIN32 |
101 | # define flukam flukam_ |
102 | # define fluka_openinp fluka_openinp_ |
103 | # define fluka_closeinp fluka_closeinp_ |
f9cb2fec |
104 | # define mcihad mcihad_ |
105 | # define mpdgha mpdgha_ |
b9d0a01d |
106 | #else |
107 | # define flukam FLUKAM |
108 | # define fluka_openinp FLUKA_OPENINP |
109 | # define fluka_closeinp FLUKA_CLOSEINP |
f9cb2fec |
110 | # define mcihad MCIHAD |
111 | # define mpdgha MPDGHA |
b9d0a01d |
112 | #endif |
113 | |
114 | extern "C" |
115 | { |
116 | // |
117 | // Prototypes for FLUKA functions |
118 | // |
119 | void type_of_call flukam(const int&); |
120 | void type_of_call fluka_openinp(const int&, DEFCHARA); |
121 | void type_of_call fluka_closeinp(const int&); |
f9cb2fec |
122 | int type_of_call mcihad(const int&); |
123 | int type_of_call mpdgha(const int&); |
b9d0a01d |
124 | } |
125 | |
126 | // |
127 | // Class implementation for ROOT |
128 | // |
03ca248b |
129 | ClassImp(TFluka) |
b9d0a01d |
130 | |
131 | // |
bf3aa28e |
132 | //---------------------------------------------------------------------------- |
133 | // TFluka constructors and destructors. |
b9d0a01d |
134 | //____________________________________________________________________________ |
135 | TFluka::TFluka() |
136 | :TVirtualMC(), |
137 | fVerbosityLevel(0), |
bf3aa28e |
138 | fInputFileName(""), |
27b2f7fe |
139 | fDetector(0), |
140 | fCurrentFlukaRegion(-1) |
b9d0a01d |
141 | { |
142 | // |
143 | // Default constructor |
144 | // |
145 | } |
146 | |
b9d0a01d |
147 | TFluka::TFluka(const char *title, Int_t verbosity) |
148 | :TVirtualMC("TFluka",title), |
149 | fVerbosityLevel(verbosity), |
bf3aa28e |
150 | fInputFileName(""), |
27b2f7fe |
151 | fDetector(0), |
152 | fCurrentFlukaRegion(-1) |
b9d0a01d |
153 | { |
154 | if (fVerbosityLevel >=3) |
155 | cout << "==> TFluka::TFluka(" << title << ") constructor called." << endl; |
156 | |
bf3aa28e |
157 | |
158 | // create geometry manager |
159 | if (fVerbosityLevel >=2) |
160 | cout << "\t* Creating G4 Geometry manager..." << endl; |
161 | fGeometryManager = new TG4GeometryManager(); |
162 | if (fVerbosityLevel >=2) |
163 | cout << "\t* Creating G4 Detector..." << endl; |
164 | fDetector = new TG4DetConstruction(); |
165 | FGeometryInit* geominit = FGeometryInit::GetInstance(); |
166 | if (geominit) |
167 | geominit->setDetConstruction(fDetector); |
168 | else { |
169 | cerr << "ERROR: Could not create FGeometryInit!" << endl; |
170 | cerr << " Exiting!!!" << endl; |
171 | abort(); |
172 | } |
173 | |
b9d0a01d |
174 | if (fVerbosityLevel >=3) |
175 | cout << "<== TFluka::TFluka(" << title << ") constructor called." << endl; |
6d4d27f2 |
176 | |
177 | fVolumeMediaMap = new TClonesArray("FlukaVolume",1000); |
178 | fNVolumes = 0; |
179 | fMediaByRegion = 0; |
b9d0a01d |
180 | } |
181 | |
bf3aa28e |
182 | TFluka::~TFluka() { |
183 | if (fVerbosityLevel >=3) |
184 | cout << "==> TFluka::~TFluka() destructor called." << endl; |
185 | |
186 | delete fGeometryManager; |
6d4d27f2 |
187 | fVolumeMediaMap->Delete(); |
188 | delete fVolumeMediaMap; |
189 | |
bf3aa28e |
190 | |
191 | if (fVerbosityLevel >=3) |
192 | cout << "<== TFluka::~TFluka() destructor called." << endl; |
193 | } |
194 | |
195 | // |
196 | //_____________________________________________________________________________ |
197 | // TFluka control methods |
b9d0a01d |
198 | //____________________________________________________________________________ |
199 | void TFluka::Init() { |
200 | if (fVerbosityLevel >=3) |
201 | cout << "==> TFluka::Init() called." << endl; |
202 | |
203 | if (fVerbosityLevel >=2) |
204 | cout << "\t* Changing lfdrtr = (" << (GLOBAL.lfdrtr?'T':'F') |
205 | << ") in fluka..." << endl; |
206 | GLOBAL.lfdrtr = true; |
207 | |
208 | if (fVerbosityLevel >=2) |
209 | cout << "\t* Opening file " << fInputFileName << endl; |
210 | const char* fname = fInputFileName; |
211 | fluka_openinp(lunin, PASSCHARA(fname)); |
212 | |
213 | if (fVerbosityLevel >=2) |
214 | cout << "\t* Calling flukam..." << endl; |
bf3aa28e |
215 | flukam(1); |
b9d0a01d |
216 | |
217 | if (fVerbosityLevel >=2) |
218 | cout << "\t* Closing file " << fInputFileName << endl; |
219 | fluka_closeinp(lunin); |
220 | |
221 | if (fVerbosityLevel >=3) |
222 | cout << "<== TFluka::Init() called." << endl; |
fa3d1cc7 |
223 | |
6d4d27f2 |
224 | FinishGeometry(); |
225 | |
b9d0a01d |
226 | } |
227 | |
bf3aa28e |
228 | void TFluka::FinishGeometry() { |
6d4d27f2 |
229 | // |
230 | // Build-up table with region to medium correspondance |
231 | // |
232 | char tmp[5]; |
233 | |
bf3aa28e |
234 | if (fVerbosityLevel >=3) |
235 | cout << "==> TFluka::FinishGeometry() called." << endl; |
236 | |
6d4d27f2 |
237 | // fGeometryManager->Ggclos(); |
bf3aa28e |
238 | |
6d4d27f2 |
239 | FGeometryInit* flugg = FGeometryInit::GetInstance(); |
240 | |
241 | fMediaByRegion = new Int_t[fNVolumes+2]; |
242 | for (Int_t i = 0; i < fNVolumes; i++) |
243 | { |
244 | FlukaVolume* vol = dynamic_cast<FlukaVolume*>((*fVolumeMediaMap)[i]); |
245 | TString volName = vol->GetName(); |
246 | Int_t media = vol->GetMedium(); |
247 | printf("Finish Geometry: volName, media %d %s %d \n", i, volName.Data(), media); |
248 | strcpy(tmp, volName.Data()); |
249 | tmp[4] = '\0'; |
b0d8df96 |
250 | flugg->SetMediumFromName(tmp, media, i+1); |
251 | fMediaByRegion[i] = media; |
27b2f7fe |
252 | } |
6d4d27f2 |
253 | |
254 | flugg->BuildMediaMap(); |
27b2f7fe |
255 | |
bf3aa28e |
256 | if (fVerbosityLevel >=3) |
257 | cout << "<== TFluka::FinishGeometry() called." << endl; |
258 | } |
259 | |
260 | void TFluka::BuildPhysics() { |
261 | if (fVerbosityLevel >=3) |
262 | cout << "==> TFluka::BuildPhysics() called." << endl; |
263 | |
264 | |
265 | if (fVerbosityLevel >=3) |
266 | cout << "<== TFluka::BuildPhysics() called." << endl; |
267 | } |
268 | |
b9d0a01d |
269 | void TFluka::ProcessEvent() { |
270 | if (fVerbosityLevel >=3) |
271 | cout << "==> TFluka::ProcessEvent() called." << endl; |
b0d8df96 |
272 | fApplication->GeneratePrimaries(); |
273 | EPISOR.lsouit = true; |
274 | flukam(1); |
b9d0a01d |
275 | if (fVerbosityLevel >=3) |
276 | cout << "<== TFluka::ProcessEvent() called." << endl; |
277 | } |
278 | |
bf3aa28e |
279 | |
b9d0a01d |
280 | void TFluka::ProcessRun(Int_t nevent) { |
281 | if (fVerbosityLevel >=3) |
282 | cout << "==> TFluka::ProcessRun(" << nevent << ") called." |
283 | << endl; |
284 | |
285 | if (fVerbosityLevel >=2) { |
286 | cout << "\t* GLOBAL.fdrtr = " << (GLOBAL.lfdrtr?'T':'F') << endl; |
287 | cout << "\t* Calling flukam again..." << endl; |
288 | } |
b0d8df96 |
289 | fApplication->InitGeometry(); |
290 | fApplication->BeginEvent(); |
291 | ProcessEvent(); |
292 | fApplication->FinishEvent(); |
b9d0a01d |
293 | if (fVerbosityLevel >=3) |
294 | cout << "<== TFluka::ProcessRun(" << nevent << ") called." |
295 | << endl; |
b0d8df96 |
296 | |
b9d0a01d |
297 | } |
298 | |
bf3aa28e |
299 | //_____________________________________________________________________________ |
300 | // methods for building/management of geometry |
301 | //____________________________________________________________________________ |
302 | // functions from GCONS |
303 | void TFluka::Gfmate(Int_t imat, char *name, Float_t &a, Float_t &z, |
304 | Float_t &dens, Float_t &radl, Float_t &absl, |
305 | Float_t* ubuf, Int_t& nbuf) { |
306 | // |
307 | fGeometryManager->Gfmate(imat, name, a, z, dens, radl, absl, ubuf, nbuf); |
308 | } |
309 | |
310 | void TFluka::Gfmate(Int_t imat, char *name, Double_t &a, Double_t &z, |
311 | Double_t &dens, Double_t &radl, Double_t &absl, |
312 | Double_t* ubuf, Int_t& nbuf) { |
313 | // |
314 | fGeometryManager->Gfmate(imat, name, a, z, dens, radl, absl, ubuf, nbuf); |
315 | } |
316 | |
317 | // detector composition |
318 | void TFluka::Material(Int_t& kmat, const char* name, Double_t a, |
319 | Double_t z, Double_t dens, Double_t radl, Double_t absl, |
320 | Float_t* buf, Int_t nwbuf) { |
321 | // |
322 | fGeometryManager |
323 | ->Material(kmat, name, a, z, dens, radl, absl, buf, nwbuf); |
324 | } |
325 | void TFluka::Material(Int_t& kmat, const char* name, Double_t a, |
326 | Double_t z, Double_t dens, Double_t radl, Double_t absl, |
327 | Double_t* buf, Int_t nwbuf) { |
328 | // |
329 | fGeometryManager |
330 | ->Material(kmat, name, a, z, dens, radl, absl, buf, nwbuf); |
331 | } |
332 | |
333 | void TFluka::Mixture(Int_t& kmat, const char *name, Float_t *a, |
334 | Float_t *z, Double_t dens, Int_t nlmat, Float_t *wmat) { |
335 | // |
336 | fGeometryManager |
337 | ->Mixture(kmat, name, a, z, dens, nlmat, wmat); |
338 | } |
339 | void TFluka::Mixture(Int_t& kmat, const char *name, Double_t *a, |
340 | Double_t *z, Double_t dens, Int_t nlmat, Double_t *wmat) { |
341 | // |
342 | fGeometryManager |
343 | ->Mixture(kmat, name, a, z, dens, nlmat, wmat); |
344 | } |
345 | |
346 | void TFluka::Medium(Int_t& kmed, const char *name, Int_t nmat, |
347 | Int_t isvol, Int_t ifield, Double_t fieldm, Double_t tmaxfd, |
348 | Double_t stemax, Double_t deemax, Double_t epsil, |
349 | Double_t stmin, Float_t* ubuf, Int_t nbuf) { |
350 | // |
351 | fGeometryManager |
352 | ->Medium(kmed, name, nmat, isvol, ifield, fieldm, tmaxfd, stemax, deemax, |
353 | epsil, stmin, ubuf, nbuf); |
354 | } |
355 | void TFluka::Medium(Int_t& kmed, const char *name, Int_t nmat, |
356 | Int_t isvol, Int_t ifield, Double_t fieldm, Double_t tmaxfd, |
357 | Double_t stemax, Double_t deemax, Double_t epsil, |
358 | Double_t stmin, Double_t* ubuf, Int_t nbuf) { |
359 | // |
360 | fGeometryManager |
361 | ->Medium(kmed, name, nmat, isvol, ifield, fieldm, tmaxfd, stemax, deemax, |
362 | epsil, stmin, ubuf, nbuf); |
363 | } |
364 | |
365 | void TFluka::Matrix(Int_t& krot, Double_t thetaX, Double_t phiX, |
366 | Double_t thetaY, Double_t phiY, Double_t thetaZ, |
367 | Double_t phiZ) { |
368 | // |
369 | fGeometryManager |
370 | ->Matrix(krot, thetaX, phiX, thetaY, phiY, thetaZ, phiZ); |
371 | } |
372 | |
373 | void TFluka::Gstpar(Int_t itmed, const char *param, Double_t parval) { |
374 | // |
375 | fGeometryManager->Gstpar(itmed, param, parval); |
376 | } |
377 | |
378 | // functions from GGEOM |
379 | Int_t TFluka::Gsvolu(const char *name, const char *shape, Int_t nmed, |
380 | Float_t *upar, Int_t np) { |
381 | // |
6d4d27f2 |
382 | // fVolumeMediaMap[TString(name)] = nmed; |
b0d8df96 |
383 | printf("TFluka::Gsvolu() name = %s, nmed = %d\n", name, nmed); |
384 | |
6d4d27f2 |
385 | TClonesArray &lvols = *fVolumeMediaMap; |
386 | new(lvols[fNVolumes++]) |
387 | FlukaVolume(name, nmed); |
388 | return fGeometryManager->Gsvolu(name, shape, nmed, upar, np); |
bf3aa28e |
389 | } |
390 | Int_t TFluka::Gsvolu(const char *name, const char *shape, Int_t nmed, |
391 | Double_t *upar, Int_t np) { |
392 | // |
6d4d27f2 |
393 | TClonesArray &lvols = *fVolumeMediaMap; |
394 | new(lvols[fNVolumes++]) |
395 | FlukaVolume(name, nmed); |
396 | |
397 | return fGeometryManager->Gsvolu(name, shape, nmed, upar, np); |
bf3aa28e |
398 | } |
399 | |
400 | void TFluka::Gsdvn(const char *name, const char *mother, Int_t ndiv, |
401 | Int_t iaxis) { |
402 | // |
b0d8df96 |
403 | // The medium of the daughter is the one of the mother |
404 | Int_t volid = TFluka::VolId(mother); |
405 | Int_t med = TFluka::VolId2Mate(volid); |
406 | TClonesArray &lvols = *fVolumeMediaMap; |
407 | new(lvols[fNVolumes++]) |
408 | FlukaVolume(name, med); |
6d4d27f2 |
409 | fGeometryManager->Gsdvn(name, mother, ndiv, iaxis); |
bf3aa28e |
410 | } |
411 | |
412 | void TFluka::Gsdvn2(const char *name, const char *mother, Int_t ndiv, |
413 | Int_t iaxis, Double_t c0i, Int_t numed) { |
414 | // |
6d4d27f2 |
415 | TClonesArray &lvols = *fVolumeMediaMap; |
416 | new(lvols[fNVolumes++]) |
417 | FlukaVolume(name, numed); |
418 | fGeometryManager->Gsdvn2(name, mother, ndiv, iaxis, c0i, numed); |
bf3aa28e |
419 | } |
420 | |
421 | void TFluka::Gsdvt(const char *name, const char *mother, Double_t step, |
422 | Int_t iaxis, Int_t numed, Int_t ndvmx) { |
6d4d27f2 |
423 | // |
424 | TClonesArray &lvols = *fVolumeMediaMap; |
425 | new(lvols[fNVolumes++]) |
426 | FlukaVolume(name, numed); |
427 | fGeometryManager->Gsdvt(name, mother, step, iaxis, numed, ndvmx); |
bf3aa28e |
428 | } |
429 | |
430 | void TFluka::Gsdvt2(const char *name, const char *mother, Double_t step, |
431 | Int_t iaxis, Double_t c0, Int_t numed, Int_t ndvmx) { |
432 | // |
6d4d27f2 |
433 | TClonesArray &lvols = *fVolumeMediaMap; |
434 | new(lvols[fNVolumes++]) |
435 | FlukaVolume(name, numed); |
436 | fGeometryManager->Gsdvt2(name, mother, step, iaxis, c0, numed, ndvmx); |
bf3aa28e |
437 | } |
438 | |
439 | void TFluka::Gsord(const char *name, Int_t iax) { |
440 | // |
441 | fGeometryManager->Gsord(name, iax); |
442 | } |
443 | |
444 | void TFluka::Gspos(const char *name, Int_t nr, const char *mother, |
445 | Double_t x, Double_t y, Double_t z, Int_t irot, |
446 | const char *konly) { |
447 | // |
448 | fGeometryManager->Gspos(name, nr, mother, x, y, z, irot, konly); |
449 | } |
450 | |
451 | void TFluka::Gsposp(const char *name, Int_t nr, const char *mother, |
452 | Double_t x, Double_t y, Double_t z, Int_t irot, |
453 | const char *konly, Float_t *upar, Int_t np) { |
454 | // |
455 | fGeometryManager->Gsposp(name, nr, mother, x, y, z, irot, konly, upar, np); |
456 | } |
457 | void TFluka::Gsposp(const char *name, Int_t nr, const char *mother, |
458 | Double_t x, Double_t y, Double_t z, Int_t irot, |
459 | const char *konly, Double_t *upar, Int_t np) { |
460 | // |
461 | fGeometryManager->Gsposp(name, nr, mother, x, y, z, irot, konly, upar, np); |
462 | } |
463 | |
464 | void TFluka::Gsbool(const char* onlyVolName, const char* manyVolName) { |
465 | // |
466 | fGeometryManager->Gsbool(onlyVolName, manyVolName); |
467 | } |
468 | |
469 | void TFluka::SetCerenkov(Int_t itmed, Int_t npckov, Float_t *ppckov, |
470 | Float_t *absco, Float_t *effic, Float_t *rindex) { |
471 | // |
472 | fGeometryManager->SetCerenkov(itmed, npckov, ppckov, absco, effic, rindex); |
473 | } |
474 | void TFluka::SetCerenkov(Int_t itmed, Int_t npckov, Double_t *ppckov, |
475 | Double_t *absco, Double_t *effic, Double_t *rindex) { |
476 | // |
477 | fGeometryManager->SetCerenkov(itmed, npckov, ppckov, absco, effic, rindex); |
478 | } |
479 | |
480 | // Euclid |
481 | void TFluka::WriteEuclid(const char* fileName, const char* topVol, |
482 | Int_t number, Int_t nlevel) { |
483 | // |
484 | fGeometryManager->WriteEuclid(fileName, topVol, number, nlevel); |
485 | } |
486 | |
487 | |
488 | |
27b2f7fe |
489 | //_____________________________________________________________________________ |
490 | // methods needed by the stepping |
491 | //____________________________________________________________________________ |
6d4d27f2 |
492 | |
27b2f7fe |
493 | Int_t TFluka::GetMedium() const { |
b0d8df96 |
494 | // |
495 | // Get the medium number for the current fluka region |
496 | // |
6d4d27f2 |
497 | FGeometryInit* flugg = FGeometryInit::GetInstance(); |
498 | return flugg->GetMedium(fCurrentFlukaRegion); |
27b2f7fe |
499 | } |
bf3aa28e |
500 | |
501 | |
502 | |
503 | //____________________________________________________________________________ |
504 | // ID <--> PDG transformations |
b9d0a01d |
505 | //_____________________________________________________________________________ |
506 | Int_t TFluka::IdFromPDG(Int_t pdg) const |
507 | { |
508 | // |
f9cb2fec |
509 | // Return Fluka code from PDG and pseudo ENDF code |
510 | |
511 | // MCIHAD() goes from pdg to fluka internal. |
512 | Int_t intfluka = mcihad(pdg); |
513 | // KPTOIP array goes from internal to official |
514 | return GetFlukaKPTOIP(intfluka); |
b9d0a01d |
515 | } |
516 | |
b9d0a01d |
517 | Int_t TFluka::PDGFromId(Int_t id) const |
518 | { |
519 | // |
f9cb2fec |
520 | // Return PDG code and pseudo ENDF code from Fluka code |
c230803a |
521 | |
bc021b12 |
522 | //IPTOKP array goes from official to internal |
b0d8df96 |
523 | if (id == 0) { |
524 | printf("PDGFromId: Error id = 0"); |
525 | return -1; |
526 | } |
527 | |
bc021b12 |
528 | Int_t intfluka = GetFlukaIPTOKP(id); |
b0d8df96 |
529 | if (intfluka == 0) { |
530 | printf("PDGFromId: Error intfluka = 0"); |
531 | return -1; |
532 | } |
533 | |
534 | //MPKDHA() goes from internal to PDG |
bc021b12 |
535 | return mpdgha(intfluka); |
6d4d27f2 |
536 | } |
537 | |
fa3d1cc7 |
538 | //_____________________________________________________________________________ |
539 | // methods for step management |
540 | //____________________________________________________________________________ |
bc021b12 |
541 | // |
542 | // set methods |
543 | // |
544 | void TFluka::SetMaxStep(Double_t) |
545 | { |
546 | // SetMaxStep is dummy procedure in TFluka ! |
547 | cout << "SetMaxStep is dummy procedure in TFluka !" << endl; |
548 | } |
549 | |
550 | void TFluka::SetMaxNStep(Int_t) |
551 | { |
552 | // SetMaxNStep is dummy procedure in TFluka ! |
553 | cout << "SetMaxNStep is dummy procedure in TFluka !" << endl; |
554 | } |
555 | |
556 | void TFluka::SetUserDecay(Int_t) |
557 | { |
558 | // SetUserDecay is dummy procedure in TFluka ! |
559 | cout << "SetUserDecay is dummy procedure in TFluka !" << endl; |
560 | } |
561 | |
fa3d1cc7 |
562 | // |
563 | // dynamic properties |
564 | // |
565 | void TFluka::TrackPosition(TLorentzVector& position) const |
566 | { |
567 | // Return the current position in the master reference frame of the |
568 | // track being transported |
569 | // TRACKR.atrack = age of the particle |
570 | // TRACKR.xtrack = x-position of the last point |
571 | // TRACKR.ytrack = y-position of the last point |
572 | // TRACKR.ztrack = z-position of the last point |
573 | position.SetX(TRACKR.xtrack[TRACKR.ntrack]); |
574 | position.SetY(TRACKR.ytrack[TRACKR.ntrack]); |
575 | position.SetZ(TRACKR.ztrack[TRACKR.ntrack]); |
576 | position.SetT(TRACKR.atrack); |
577 | } |
578 | |
579 | void TFluka::TrackMomentum(TLorentzVector& momentum) const |
580 | { |
581 | // Return the direction and the momentum (GeV/c) of the track |
582 | // currently being transported |
583 | // TRACKR.ptrack = momentum of the particle (not always defined, if |
584 | // < 0 must be obtained from etrack) |
585 | // TRACKR.cx,y,ztrck = direction cosines of the current particle |
586 | // TRACKR.etrack = total energy of the particle |
587 | // TRACKR.jtrack = identity number of the particle |
588 | // PAPROP.am[TRACKR.jtrack] = particle mass in gev |
589 | if (TRACKR.ptrack >= 0) { |
590 | momentum.SetPx(TRACKR.ptrack*TRACKR.cxtrck); |
591 | momentum.SetPy(TRACKR.ptrack*TRACKR.cytrck); |
592 | momentum.SetPz(TRACKR.ptrack*TRACKR.cztrck); |
593 | momentum.SetE(TRACKR.etrack); |
594 | return; |
595 | } |
596 | else { |
b8b430a9 |
597 | Double_t p = sqrt(TRACKR.etrack*TRACKR.etrack - PAPROP.am[TRACKR.jtrack+6]*PAPROP.am[TRACKR.jtrack+6]); |
fa3d1cc7 |
598 | momentum.SetPx(p*TRACKR.cxtrck); |
599 | momentum.SetPy(p*TRACKR.cytrck); |
600 | momentum.SetPz(p*TRACKR.cztrck); |
601 | momentum.SetE(TRACKR.etrack); |
602 | return; |
603 | } |
604 | } |
605 | |
606 | Double_t TFluka::TrackStep() const |
607 | { |
608 | // Return the length in centimeters of the current step |
609 | // TRACKR.ctrack = total curved path |
610 | return TRACKR.ctrack; |
611 | } |
612 | |
613 | Double_t TFluka::TrackLength() const |
614 | { |
b8b430a9 |
615 | // Still wrong !!! |
616 | // This is the sum of substeps !!! |
617 | // TRACKR.ctrack = total curved path of the current step |
618 | // Sum of the substeps is identical to TRACKR.ctrack if the is no mag. field |
619 | // The sum of all step length starting from the beginning of the track |
fa3d1cc7 |
620 | // for the time being returns only the length in centimeters of the current step |
b8b430a9 |
621 | Double_t sum = 0; |
622 | for ( Int_t j=0;j<TRACKR.ntrack;j++) { |
623 | sum +=TRACKR.ttrack[j]; |
624 | } |
625 | return sum; |
fa3d1cc7 |
626 | } |
627 | |
628 | Double_t TFluka::TrackTime() const |
629 | { |
630 | // Return the current time of flight of the track being transported |
631 | // TRACKR.atrack = age of the particle |
632 | return TRACKR.atrack; |
633 | } |
634 | |
635 | Double_t TFluka::Edep() const |
636 | { |
637 | // Energy deposition |
638 | // if TRACKR.ntrack = 0, TRACKR.mtrack = 0: |
639 | // -->local energy deposition (the value and the point are not recorded in TRACKR) |
640 | // but in the variable "rull" of the procedure "endraw.cxx" |
641 | // if TRACKR.ntrack > 0, TRACKR.mtrack = 0: |
642 | // -->no energy loss along the track |
643 | // if TRACKR.ntrack > 0, TRACKR.mtrack > 0: |
644 | // -->energy loss distributed along the track |
645 | // TRACKR.dtrack = energy deposition of the jth deposition even |
646 | if (TRACKR.ntrack == 0 && TRACKR.mtrack == 0) |
647 | return fRull; |
648 | else { |
649 | Double_t sum = 0; |
650 | for ( Int_t j=0;j<TRACKR.mtrack;j++) { |
651 | sum +=TRACKR.dtrack[j]; |
652 | } |
653 | return sum; |
654 | } |
655 | } |
656 | |
657 | Int_t TFluka::TrackPid() const |
658 | { |
659 | // Return the id of the particle transported |
660 | // TRACKR.jtrack = identity number of the particle |
661 | return PDGFromId(TRACKR.jtrack); |
662 | } |
663 | |
664 | Double_t TFluka::TrackCharge() const |
665 | { |
666 | // Return charge of the track currently transported |
667 | // PAPROP.ichrge = electric charge of the particle |
bc021b12 |
668 | // TRACKR.jtrack = identity number of the particle |
fa3d1cc7 |
669 | return PAPROP.ichrge[TRACKR.jtrack+6]; |
670 | } |
671 | |
672 | Double_t TFluka::TrackMass() const |
673 | { |
674 | // PAPROP.am = particle mass in GeV |
bc021b12 |
675 | // TRACKR.jtrack = identity number of the particle |
fa3d1cc7 |
676 | return PAPROP.am[TRACKR.jtrack+6]; |
677 | } |
678 | |
679 | Double_t TFluka::Etot() const |
680 | { |
681 | // TRACKR.etrack = total energy of the particle |
682 | return TRACKR.etrack; |
683 | } |
684 | |
685 | // |
686 | // track status |
687 | // |
688 | Bool_t TFluka::IsNewTrack() const |
689 | { |
690 | // ???????????????, |
691 | // True if the track is not at the boundary of the current volume |
692 | // Not true in some cases in bxdraw - to be solved |
693 | return 1; |
694 | } |
695 | |
696 | Bool_t TFluka::IsTrackInside() const |
697 | { |
698 | // True if the track is not at the boundary of the current volume |
699 | // In Fluka a step is always inside one kind of material |
700 | // If the step would go behind the region of one material, |
701 | // it will be shortened to reach only the boundary. |
702 | // Therefore IsTrackInside() is always true. |
703 | // Not true in some cases in bxdraw - to be solved |
704 | return 1; |
705 | } |
706 | |
707 | Bool_t TFluka::IsTrackEntering() const |
708 | { |
709 | // True if this is the first step of the track in the current volume |
710 | // Boundary- (X) crossing |
711 | // Icode = 19: boundary crossing - call from Kaskad |
712 | // Icode = 29: boundary crossing - call from Emfsco |
713 | // Icode = 39: boundary crossing - call from Kasneu |
714 | // Icode = 49: boundary crossing - call from Kashea |
715 | // Icode = 59: boundary crossing - call from Kasoph |
716 | if (fIcode == 19 || |
717 | fIcode == 29 || |
718 | fIcode == 39 || |
719 | fIcode == 49 || |
720 | fIcode == 59) return 1; |
721 | else return 0; |
722 | } |
723 | |
724 | Bool_t TFluka::IsTrackExiting() const |
725 | { |
726 | // True if this is the last step of the track in the current volume |
727 | // Boundary- (X) crossing |
728 | // Icode = 19: boundary crossing - call from Kaskad |
729 | // Icode = 29: boundary crossing - call from Emfsco |
730 | // Icode = 39: boundary crossing - call from Kasneu |
731 | // Icode = 49: boundary crossing - call from Kashea |
732 | // Icode = 59: boundary crossing - call from Kasoph |
733 | if (fIcode == 19 || |
734 | fIcode == 29 || |
735 | fIcode == 39 || |
736 | fIcode == 49 || |
737 | fIcode == 59) return 1; |
738 | else return 0; |
739 | } |
740 | |
741 | Bool_t TFluka::IsTrackOut() const |
742 | { |
743 | // True if the track is out of the setup |
744 | // means escape |
745 | // Icode = 14: escape - call from Kaskad |
746 | // Icode = 23: escape - call from Emfsco |
747 | // Icode = 32: escape - call from Kasneu |
748 | // Icode = 40: escape - call from Kashea |
749 | // Icode = 51: escape - call from Kasoph |
750 | if (fIcode == 14 || |
751 | fIcode == 23 || |
752 | fIcode == 32 || |
753 | fIcode == 40 || |
754 | fIcode == 51) return 1; |
755 | else return 0; |
756 | } |
757 | |
758 | Bool_t TFluka::IsTrackDisappeared() const |
759 | { |
760 | // means all inelastic interactions and decays |
761 | // fIcode from usdraw |
762 | if (fIcode == 101 || // inelastic interaction |
763 | fIcode == 102 || // particle decay |
764 | fIcode == 214 || // in-flight annihilation |
765 | fIcode == 215 || // annihilation at rest |
766 | fIcode == 217 || // pair production |
767 | fIcode == 221) return 1; |
768 | else return 0; |
769 | } |
770 | |
771 | Bool_t TFluka::IsTrackStop() const |
772 | { |
773 | // True if the track energy has fallen below the threshold |
774 | // means stopped by signal or below energy threshold |
775 | // Icode = 12: stopping particle - call from Kaskad |
776 | // Icode = 15: time kill - call from Kaskad |
777 | // Icode = 21: below threshold, iarg=1 - call from Emfsco |
778 | // Icode = 22: below threshold, iarg=2 - call from Emfsco |
779 | // Icode = 24: time kill - call from Emfsco |
780 | // Icode = 31: below threshold - call from Kasneu |
781 | // Icode = 33: time kill - call from Kasneu |
782 | // Icode = 41: time kill - call from Kashea |
783 | // Icode = 52: time kill - call from Kasoph |
784 | if (fIcode == 12 || |
785 | fIcode == 15 || |
786 | fIcode == 21 || |
787 | fIcode == 22 || |
788 | fIcode == 24 || |
789 | fIcode == 31 || |
790 | fIcode == 33 || |
791 | fIcode == 41 || |
792 | fIcode == 52) return 1; |
793 | else return 0; |
794 | } |
795 | |
796 | Bool_t TFluka::IsTrackAlive() const |
797 | { |
798 | // means not disappeared or not out |
799 | if (IsTrackDisappeared() || IsTrackOut() ) return 0; |
800 | else return 1; |
801 | } |
802 | |
803 | // |
804 | // secondaries |
805 | // |
806 | |
807 | Int_t TFluka::NSecondaries() const |
808 | // Number of secondary particles generated in the current step |
bc021b12 |
809 | // FINUC.np = number of secondaries except light and heavy ions |
b8b430a9 |
810 | // FHEAVY.npheav = number of secondaries for light and heavy secondary ions |
fa3d1cc7 |
811 | { |
bc021b12 |
812 | return FINUC.np + FHEAVY.npheav; |
fa3d1cc7 |
813 | } |
814 | |
815 | void TFluka::GetSecondary(Int_t isec, Int_t& particleId, |
816 | TLorentzVector& position, TLorentzVector& momentum) |
fa3d1cc7 |
817 | { |
bc021b12 |
818 | if (isec >= 0 && isec < FINUC.np) { |
b8b430a9 |
819 | // more fine condition depending on icode |
820 | // icode = 100 ? |
821 | // icode = 101 OK |
822 | // icode = 102 OK |
823 | // icode = 103 ? |
824 | // icode = 104 ? |
825 | // icode = 105 ? |
826 | // icode = 208 ? |
827 | // icode = 210 ? |
828 | // icode = 212 ? |
829 | // icode = 214 OK |
830 | // icode = 215 OK |
831 | // icode = 219 ? |
832 | // icode = 221 OK |
833 | // icode = 225 ? |
834 | // icode = 300 ? |
835 | // icode = 400 ? |
836 | |
bc021b12 |
837 | particleId = PDGFromId(FINUC.kpart[isec]); |
838 | position.SetX(fXsco); |
839 | position.SetY(fYsco); |
840 | position.SetZ(fZsco); |
841 | position.SetT(TRACKR.atrack); |
842 | // position.SetT(TRACKR.atrack+FINUC.agesec[isec]); //not yet implem. |
843 | momentum.SetPx(FINUC.plr[isec]*FINUC.cxr[isec]); |
844 | momentum.SetPy(FINUC.plr[isec]*FINUC.cyr[isec]); |
845 | momentum.SetPz(FINUC.plr[isec]*FINUC.czr[isec]); |
846 | momentum.SetE(FINUC.tki[isec] + PAPROP.am[FINUC.kpart[isec]+6]); |
fa3d1cc7 |
847 | } |
bc021b12 |
848 | if (isec >= FINUC.np && isec < FINUC.np + FHEAVY.npheav) { |
849 | Int_t jsec = isec - FINUC.np; |
850 | particleId = FHEAVY.kheavy[jsec]; // this is Fluka id !!! |
851 | position.SetX(fXsco); |
852 | position.SetY(fYsco); |
853 | position.SetZ(fZsco); |
854 | position.SetT(TRACKR.atrack); |
855 | // position.SetT(TRACKR.atrack+FHEAVY.agheav[jsec]); //not yet implem. |
856 | momentum.SetPx(FHEAVY.pheavy[jsec]*FHEAVY.cxheav[jsec]); |
857 | momentum.SetPy(FHEAVY.pheavy[jsec]*FHEAVY.cyheav[jsec]); |
858 | momentum.SetPz(FHEAVY.pheavy[jsec]*FHEAVY.czheav[jsec]); |
859 | if (FHEAVY.tkheav[jsec] >= 3 && FHEAVY.tkheav[jsec] <= 6) |
860 | momentum.SetE(FHEAVY.tkheav[jsec] + PAPROP.am[jsec+6]); |
861 | else if (FHEAVY.tkheav[jsec] > 6) |
862 | momentum.SetE(FHEAVY.tkheav[jsec] + FHEAVY.amnhea[jsec]); // to be checked !!! |
863 | } |
fa3d1cc7 |
864 | } |
865 | |
bc021b12 |
866 | TMCProcess TFluka::ProdProcess(Int_t isec) const |
fa3d1cc7 |
867 | // Name of the process that has produced the secondary particles |
868 | // in the current step |
bc021b12 |
869 | { |
870 | const TMCProcess kIpNoProc = kPNoProcess; |
871 | const TMCProcess kIpPDecay = kPDecay; |
872 | const TMCProcess kIpPPair = kPPair; |
873 | //const TMCProcess kIpPPairFromPhoton = kPPairFromPhoton; |
874 | //const TMCProcess kIpPPairFromVirtualPhoton = kPPairFromVirtualPhoton; |
875 | const TMCProcess kIpPCompton = kPCompton; |
876 | const TMCProcess kIpPPhotoelectric = kPPhotoelectric; |
877 | const TMCProcess kIpPBrem = kPBrem; |
878 | //const TMCProcess kIpPBremFromHeavy = kPBremFromHeavy; |
879 | //const TMCProcess kIpPBremFromElectronOrPositron = kPBremFromElectronOrPositron; |
880 | const TMCProcess kIpPDeltaRay = kPDeltaRay; |
881 | //const TMCProcess kIpPMoller = kPMoller; |
882 | //const TMCProcess kIpPBhabha = kPBhabha; |
883 | const TMCProcess kIpPAnnihilation = kPAnnihilation; |
884 | //const TMCProcess kIpPAnnihilInFlight = kPAnnihilInFlight; |
885 | //const TMCProcess kIpPAnnihilAtRest = kPAnnihilAtRest; |
886 | const TMCProcess kIpPHadronic = kPHadronic; |
887 | const TMCProcess kIpPMuonNuclear = kPMuonNuclear; |
888 | const TMCProcess kIpPPhotoFission = kPPhotoFission; |
889 | const TMCProcess kIpPRayleigh = kPRayleigh; |
b0d8df96 |
890 | // const TMCProcess kIpPCerenkov = kPCerenkov; |
891 | // const TMCProcess kIpPSynchrotron = kPSynchrotron; |
bc021b12 |
892 | |
893 | Int_t mugamma = TRACKR.jtrack == 7 || TRACKR.jtrack == 10 || TRACKR.jtrack == 11; |
894 | if (fIcode == 102) return kIpPDecay; |
895 | else if (fIcode == 104 || fIcode == 217) return kIpPPair; |
896 | //else if (fIcode == 104) return kIpPairFromPhoton; |
897 | //else if (fIcode == 217) return kIpPPairFromVirtualPhoton; |
898 | else if (fIcode == 219) return kIpPCompton; |
899 | else if (fIcode == 221) return kIpPPhotoelectric; |
900 | else if (fIcode == 105 || fIcode == 208) return kIpPBrem; |
901 | //else if (fIcode == 105) return kIpPBremFromHeavy; |
902 | //else if (fIcode == 208) return kPBremFromElectronOrPositron; |
903 | else if (fIcode == 103 || fIcode == 400) return kIpPDeltaRay; |
904 | else if (fIcode == 210 || fIcode == 212) return kIpPDeltaRay; |
905 | //else if (fIcode == 210) return kIpPMoller; |
906 | //else if (fIcode == 212) return kIpPBhabha; |
907 | else if (fIcode == 214 || fIcode == 215) return kIpPAnnihilation; |
908 | //else if (fIcode == 214) return kIpPAnnihilInFlight; |
909 | //else if (fIcode == 215) return kIpPAnnihilAtRest; |
910 | else if (fIcode == 101) return kIpPHadronic; |
911 | else if (fIcode == 101) { |
912 | if (!mugamma) return kIpPHadronic; |
913 | else if (TRACKR.jtrack == 7) return kIpPPhotoFission; |
914 | else return kIpPMuonNuclear; |
915 | } |
916 | else if (fIcode == 225) return kIpPRayleigh; |
917 | // Fluka codes 100, 300 and 400 still to be investigasted |
918 | else return kIpNoProc; |
919 | } |
fa3d1cc7 |
920 | |
921 | //Int_t StepProcesses(TArrayI &proc) const |
922 | // Return processes active in the current step |
923 | //{ |
924 | //ck = total energy of the particl ???????????????? |
925 | //} |
926 | |
927 | |
b0d8df96 |
928 | Int_t TFluka::VolId2Mate(Int_t id) const |
929 | { |
930 | // |
931 | // Returns the material number for a given volume ID |
932 | // |
933 | printf("VolId2Mate %d %d\n", id, fMediaByRegion[id]); |
934 | return fMediaByRegion[id-1]; |
935 | } |
936 | |
937 | const char* TFluka::VolName(Int_t id) const |
938 | { |
939 | // |
940 | // Returns the volume name for a given volume ID |
941 | // |
942 | FlukaVolume* vol = dynamic_cast<FlukaVolume*>((*fVolumeMediaMap)[id-1]); |
943 | const char* name = vol->GetName(); |
944 | printf("VolName %d %s \n", id, name); |
945 | return name; |
946 | } |
947 | |
948 | Int_t TFluka::VolId(const Text_t* volName) const |
949 | { |
950 | // |
951 | // Converts from volume name to volume ID. |
952 | // Time consuming. (Only used during set-up) |
953 | // Could be replaced by hash-table |
954 | // |
955 | char tmp[5]; |
956 | Int_t i =0; |
957 | for (i = 0; i < fNVolumes; i++) |
958 | { |
959 | FlukaVolume* vol = dynamic_cast<FlukaVolume*>((*fVolumeMediaMap)[i]); |
960 | TString name = vol->GetName(); |
961 | strcpy(tmp, name.Data()); |
962 | tmp[4] = '\0'; |
963 | if (!strcmp(tmp, volName)) break; |
964 | } |
965 | i++; |
966 | |
967 | return i; |
968 | } |
969 | |
970 | |
971 | Int_t TFluka::CurrentVolID(Int_t& copyNo) const |
972 | { |
973 | // |
974 | // Return the logical id and copy number corresponding to the current fluka region |
975 | // |
976 | int ir = fCurrentFlukaRegion; |
977 | int id = (FGeometryInit::GetInstance())->CurrentVolID(ir, copyNo); |
978 | printf("CurrentVolID: %d %d %d \n", ir, id, copyNo); |
979 | return id; |
980 | |
981 | } |
982 | |
983 | Int_t TFluka::CurrentVolOffID(Int_t off, Int_t& copyNo) const |
984 | { |
985 | // |
986 | // Return the logical id and copy number of off'th mother |
987 | // corresponding to the current fluka region |
988 | // |
989 | if (off == 0) |
990 | return CurrentVolID(copyNo); |
991 | |
992 | int ir = fCurrentFlukaRegion; |
993 | int id = (FGeometryInit::GetInstance())->CurrentVolOffID(ir, off, copyNo); |
994 | |
995 | printf("CurrentVolOffID: %d %d %d \n", ir, id, copyNo); |
996 | if (id == -1) |
997 | printf("CurrentVolOffID: Warning Mother not found !!!\n"); |
998 | return id; |
999 | } |
1000 | |
1001 | |
1002 | const char* TFluka::CurrentVolName() const |
1003 | { |
1004 | // |
1005 | // Return the current volume name |
1006 | // |
1007 | Int_t copy; |
1008 | Int_t id = TFluka::CurrentVolID(copy); |
1009 | const char* name = TFluka::VolName(id); |
1010 | printf("CurrentVolumeName: %d %s \n", fCurrentFlukaRegion, name); |
1011 | return name; |
1012 | } |
1013 | |
1014 | const char* TFluka::CurrentVolOffName(Int_t off) const |
1015 | { |
1016 | // |
1017 | // Return the volume name of the off'th mother of the current volume |
1018 | // |
1019 | Int_t copy; |
1020 | Int_t id = TFluka::CurrentVolOffID(off, copy); |
1021 | const char* name = TFluka::VolName(id); |
1022 | printf("CurrentVolumeOffName: %d %s \n", fCurrentFlukaRegion, name); |
1023 | return name; |
1024 | } |
1025 | |
1026 | Int_t TFluka::CurrentMaterial(Float_t &a, Float_t &z, |
1027 | Float_t &dens, Float_t &radl, Float_t &absl) const |
1028 | { |
1029 | // |
1030 | // Return the current medium number |
1031 | // |
1032 | Int_t copy; |
1033 | Int_t id = TFluka::CurrentVolID(copy); |
1034 | Int_t med = TFluka::VolId2Mate(id); |
1035 | printf("CurrentMaterial: %d %d \n", fCurrentFlukaRegion, med); |
1036 | return med; |
1037 | } |
1038 | |
1039 | |
fa3d1cc7 |
1040 | // =============================================================== |
1041 | void TFluka::FutoTest() |
1042 | { |
1043 | Int_t icode, mreg, newreg, particleId; |
1044 | // Int_t medium; |
1045 | Double_t rull, xsco, ysco, zsco; |
1046 | TLorentzVector position, momentum; |
1047 | icode = GetIcode(); |
1048 | if (icode == 0) { |
1049 | cout << " icode=" << icode << endl; |
1050 | /* |
1051 | cout << "TLorentzVector positionX=" << position.X() |
1052 | << "positionY=" << position.Y() |
1053 | << "positionZ=" << position.Z() |
1054 | << "timeT=" << position.T() << endl; |
1055 | cout << "TLorentzVector momentumX=" << momentum.X() |
1056 | << "momentumY=" << momentum.Y() |
1057 | << "momentumZ=" << momentum.Z() |
1058 | << "energyE=" << momentum.E() << endl; |
1059 | cout << "TrackPid=" << TrackPid() << endl; |
1060 | */ |
1061 | } |
1062 | |
1063 | else if (icode > 0 && icode <= 5) { |
1064 | // mgdraw |
1065 | mreg = GetMreg(); |
1066 | // medium = GetMedium(); |
1067 | cout << " icode=" << icode |
1068 | << " mreg=" << mreg |
1069 | // << " medium=" << medium |
1070 | << endl; |
1071 | TrackPosition(position); |
1072 | TrackMomentum(momentum); |
1073 | cout << "TLorentzVector positionX=" << position.X() |
1074 | << "positionY=" << position.Y() |
1075 | << "positionZ=" << position.Z() |
1076 | << "timeT=" << position.T() << endl; |
1077 | cout << "TLorentzVector momentumX=" << momentum.X() |
1078 | << "momentumY=" << momentum.Y() |
1079 | << "momentumZ=" << momentum.Z() |
1080 | << "energyE=" << momentum.E() << endl; |
1081 | cout << "TrackStep=" << TrackStep() << endl; |
1082 | cout << "TrackLength=" << TrackLength() << endl; |
1083 | cout << "TrackTime=" << TrackTime() << endl; |
1084 | cout << "Edep=" << Edep() << endl; |
1085 | cout << "TrackPid=" << TrackPid() << endl; |
1086 | cout << "TrackCharge=" << TrackCharge() << endl; |
1087 | cout << "TrackMass=" << TrackMass() << endl; |
1088 | cout << "Etot=" << Etot() << endl; |
1089 | cout << "IsNewTrack=" << IsNewTrack() << endl; |
1090 | cout << "IsTrackInside=" << IsTrackInside() << endl; |
1091 | cout << "IsTrackEntering=" << IsTrackEntering() << endl; |
1092 | cout << "IsTrackExiting=" << IsTrackExiting() << endl; |
1093 | cout << "IsTrackOut=" << IsTrackOut() << endl; |
1094 | cout << "IsTrackDisappeared=" << IsTrackDisappeared() << endl; |
1095 | cout << "IsTrackAlive=" << IsTrackAlive() << endl; |
1096 | } |
1097 | |
1098 | else if((icode >= 10 && icode <= 15) || |
1099 | (icode >= 20 && icode <= 24) || |
1100 | (icode >= 30 && icode <= 33) || |
1101 | (icode >= 40 && icode <= 41) || |
1102 | (icode >= 50 && icode <= 52)) { |
1103 | // endraw |
1104 | mreg = GetMreg(); |
1105 | // medium = GetMedium(); |
1106 | rull = GetRull(); |
1107 | xsco = GetXsco(); |
1108 | ysco = GetYsco(); |
1109 | zsco = GetZsco(); |
1110 | cout << " icode=" << icode |
1111 | << " mreg=" << mreg |
1112 | // << " medium=" << medium |
1113 | << " rull=" << rull |
1114 | << " xsco=" << xsco |
1115 | << " ysco=" << ysco |
1116 | << " zsco=" << zsco << endl; |
1117 | TrackPosition(position); |
1118 | TrackMomentum(momentum); |
1119 | cout << "Edep=" << Edep() << endl; |
1120 | cout << "Etot=" << Etot() << endl; |
1121 | cout << "TrackPid=" << TrackPid() << endl; |
1122 | cout << "TrackCharge=" << TrackCharge() << endl; |
1123 | cout << "TrackMass=" << TrackMass() << endl; |
1124 | cout << "IsTrackOut=" << IsTrackOut() << endl; |
1125 | cout << "IsTrackDisappeared=" << IsTrackDisappeared() << endl; |
1126 | cout << "IsTrackStop=" << IsTrackStop() << endl; |
1127 | cout << "IsTrackAlive=" << IsTrackAlive() << endl; |
1128 | } |
1129 | |
1130 | else if((icode >= 100 && icode <= 105) || |
1131 | (icode == 208) || |
1132 | (icode == 210) || |
1133 | (icode == 212) || |
1134 | (icode >= 214 && icode <= 215) || |
1135 | (icode == 217) || |
1136 | (icode == 219) || |
1137 | (icode == 221) || |
1138 | (icode == 225) || |
1139 | (icode == 300) || |
1140 | (icode == 400)) { |
1141 | // usdraw |
1142 | mreg = GetMreg(); |
1143 | // medium = GetMedium(); |
1144 | xsco = GetXsco(); |
1145 | ysco = GetYsco(); |
1146 | zsco = GetZsco(); |
1147 | cout << " icode=" << icode |
1148 | << " mreg=" << mreg |
1149 | // << " medium=" << medium |
1150 | << " xsco=" << xsco |
1151 | << " ysco=" << ysco |
1152 | << " zsco=" << zsco << endl; |
1153 | cout << "TrackPid=" << TrackPid() << endl; |
1154 | cout << "NSecondaries=" << NSecondaries() << endl; |
1155 | for (Int_t isec=0; isec< NSecondaries(); isec++) { |
fa3d1cc7 |
1156 | TFluka::GetSecondary(isec, particleId, position, momentum); |
1157 | cout << "TLorentzVector positionX=" << position.X() |
1158 | << "positionY=" << position.Y() |
1159 | << "positionZ=" << position.Z() |
1160 | << "timeT=" << position.T() << endl; |
1161 | cout << "TLorentzVector momentumX=" << momentum.X() |
1162 | << "momentumY=" << momentum.Y() |
1163 | << "momentumZ=" << momentum.Z() |
1164 | << "energyE=" << momentum.E() << endl; |
1165 | cout << "TrackPid=" << particleId << endl; |
1166 | |
1167 | } |
1168 | } |
1169 | |
1170 | else if((icode == 19) || |
1171 | (icode == 29) || |
1172 | (icode == 39) || |
1173 | (icode == 49) || |
1174 | (icode == 59)) { |
1175 | mreg = GetMreg(); |
1176 | // medium = GetMedium(); |
1177 | newreg = GetNewreg(); |
1178 | xsco = GetXsco(); |
1179 | ysco = GetYsco(); |
1180 | zsco = GetZsco(); |
1181 | cout << " icode=" << icode |
1182 | << " mreg=" << mreg |
1183 | // << " medium=" << medium |
1184 | << " newreg=" << newreg |
1185 | << " xsco=" << xsco |
1186 | << " ysco=" << ysco |
1187 | << " zsco=" << zsco << endl; |
1188 | } |
1189 | // |
1190 | // ==================================================================== |
1191 | // |
1192 | |
1193 | |
1194 | |
1195 | } // end of FutoTest |
b0d8df96 |
1196 | |