]>
Commit | Line | Data |
---|---|---|
1 | #include <Riostream.h> | |
2 | #include "TVirtualMCApplication.h" | |
3 | #include "TVirtualMCStack.h" | |
4 | #include "TFluka.h" | |
5 | #include "TFlukaCodes.h" | |
6 | ||
7 | // Fluka includes | |
8 | #include "Fdimpar.h" //(DIMPAR) | |
9 | #include "Fdblprc.h" //(DBLPRC) | |
10 | #include "Ftrackr.h" //(TRACKR) | |
11 | #include "Fopphst.h" //(OPPHST) | |
12 | #include "Fflkstk.h" //(FLKSTK) | |
13 | #include "Fltclcm.h" //(LTCLCM) | |
14 | #include "Fpaprop.h" //(PAPROP) | |
15 | #include "Falldlt.h" //(ALLDLT) | |
16 | #include "Fdpdxcm.h" //(DPDXCM) | |
17 | #include "Fflkmat.h" //(FLKMAT) | |
18 | ||
19 | #include "TGeoManager.h" | |
20 | ||
21 | #ifndef WIN32 | |
22 | # define mgdraw mgdraw_ | |
23 | #else | |
24 | # define mgdraw MGDRAW | |
25 | #endif | |
26 | ||
27 | extern "C" { | |
28 | void mgdraw(Int_t& icode, Int_t& mreg) | |
29 | { | |
30 | TFluka* fluka = (TFluka*) gMC; | |
31 | if (mreg == fluka->GetDummyRegion()) return; | |
32 | // | |
33 | // Make sure that stack has currrent track Id | |
34 | // | |
35 | Int_t trackId = -1; | |
36 | TVirtualMCStack* cppstack = fluka->GetStack(); | |
37 | ||
38 | if (TRACKR.jtrack == -1) { | |
39 | trackId = OPPHST.louopp[OPPHST.lstopp]; | |
40 | if (trackId == 0) { | |
41 | trackId = FLKSTK.ispark[FLKSTK.npflka][mkbmx2-1]; | |
42 | } | |
43 | } else { | |
44 | trackId = TRACKR.ispusr[mkbmx2-1]; | |
45 | } | |
46 | ||
47 | Int_t verbosityLevel = fluka->GetVerbosityLevel(); | |
48 | if (TRACKR.jtrack < -6) { | |
49 | // from -7 to -12 = "heavy" fragment | |
50 | // assing parent id | |
51 | // id < -6 was skipped in stuprf => if (kpart < -6) return; | |
52 | if (verbosityLevel >= 3) { | |
53 | cout << "mgdraw: (heavy fragment) jtrack < -6 =" << TRACKR.jtrack | |
54 | << " assign parent pdg=" << fluka->PDGFromId(TRACKR.ispusr[mkbmx2 - 3]) << endl; | |
55 | } | |
56 | // TRACKR.jtrack = TRACKR.ispusr[mkbmx2 - 3]; | |
57 | } | |
58 | ||
59 | cppstack->SetCurrentTrack(trackId); | |
60 | // | |
61 | // | |
62 | Int_t mlttc = TRACKR.lt1trk; // LTCLCM.mlatm1; | |
63 | fluka->SetMreg(mreg, mlttc); | |
64 | fluka->SetIcode((FlukaProcessCode_t) icode); | |
65 | fluka->SetCaller(kMGDRAW); | |
66 | ||
67 | Int_t nodeId; | |
68 | Int_t volId = fluka->CurrentVolID(nodeId); | |
69 | Int_t crtlttc = gGeoManager->GetCurrentNodeId()+1; | |
70 | ||
71 | // check region lattice consistency (debug Ernesto) | |
72 | // ***************************************************** | |
73 | if(verbosityLevel>=3 && mreg != volId && !gGeoManager->IsOutside() ) { | |
74 | cout << " mgdraw: track=" << trackId << " pdg=" << fluka->PDGFromId(TRACKR.jtrack) | |
75 | << " icode=" << icode << " gNstep=" << fluka->GetNstep() << endl | |
76 | << " fluka mreg=" << mreg << " mlttc=" << mlttc << endl | |
77 | << " TGeo volId=" << volId << " crtlttc=" << crtlttc << endl | |
78 | << " common TRACKR lt1trk=" << TRACKR.lt1trk << " lt2trk=" << TRACKR.lt2trk << endl | |
79 | << " common LTCLCM newlat=" << LTCLCM.newlat << " mlatld=" << LTCLCM.mlatld << endl | |
80 | << " mlatm1=" << LTCLCM.mlatm1 << " mltsen=" << LTCLCM.mltsen << endl | |
81 | << " mltsm1=" << LTCLCM.mltsm1 << " mlattc=" << LTCLCM.mlattc << endl; | |
82 | if( mlttc == crtlttc ) cout << " *************************************************************" << endl; | |
83 | } | |
84 | // ***************************************************** | |
85 | ||
86 | ||
87 | Int_t med = FLKMAT.medium[mreg - 1]; // Medium | |
88 | Int_t msd = DPDXCM.msdpdx[med - 1]; // Iionisation model | |
89 | ||
90 | if (!TRACKR.ispusr[mkbmx2 - 2]) { | |
91 | if (verbosityLevel >= 3) { | |
92 | cout << endl << "mgdraw: energy deposition for:" << trackId | |
93 | << " icode=" << icode | |
94 | << " pdg=" << fluka->PDGFromId(TRACKR.jtrack) | |
95 | << " flukaid="<< TRACKR.jtrack | |
96 | << " mreg=" << mreg | |
97 | << " np =" << ALLDLT.nalldl | |
98 | << endl; | |
99 | ||
100 | } | |
101 | ||
102 | if (msd > 0) { | |
103 | // | |
104 | // Primary ionsiations | |
105 | Int_t nprim = ALLDLT.nalldl; | |
106 | if (nprim >= mxalld) { | |
107 | nprim = mxalld; | |
108 | Warning("mgdraw", "nprim > mxalld, nprim: %6d pdg: %6d mreg %6d p %13.3f step %13.3f\n", | |
109 | ALLDLT.nalldl, | |
110 | fluka->PDGFromId(TRACKR.jtrack), | |
111 | mreg, | |
112 | TRACKR.ptrack, | |
113 | TRACKR.ctrack); | |
114 | } | |
115 | fluka->PrimaryIonisationStepping(nprim); | |
116 | } else { | |
117 | // Single step | |
118 | (TVirtualMCApplication::Instance())->Stepping(); | |
119 | fluka->SetTrackIsNew(kFALSE); | |
120 | } | |
121 | ||
122 | } else { | |
123 | // | |
124 | // Tracking is being resumed after secondary tracking | |
125 | // | |
126 | if (verbosityLevel >= 3) { | |
127 | cout << endl << "mgdraw: resuming Stepping(): " << trackId << endl; | |
128 | } | |
129 | ||
130 | fluka->SetTrackIsNew(kTRUE); | |
131 | fluka->SetCaller(kMGResumedTrack); | |
132 | (TVirtualMCApplication::Instance())->Stepping(); | |
133 | ||
134 | // Reset flag and stored values | |
135 | TRACKR.ispusr[mkbmx2 - 2] = 0; | |
136 | for (Int_t i = 0; i < 9; i++) TRACKR.spausr[i] = -1.; | |
137 | ||
138 | ||
139 | if (verbosityLevel >= 3) { | |
140 | cout << endl << " !!! I am in mgdraw - first Stepping() after resume: " << icode << endl; | |
141 | cout << " Track= " << trackId << " region = " << mreg << endl; | |
142 | } | |
143 | ||
144 | fluka->SetTrackIsNew(kFALSE); | |
145 | fluka->SetCaller(kMGDRAW); | |
146 | if (msd == 0) { | |
147 | (TVirtualMCApplication::Instance())->Stepping(); | |
148 | } else { | |
149 | Int_t nprim = ALLDLT.nalldl; | |
150 | // Protection against nprim > mxalld | |
151 | if (nprim >= mxalld) { | |
152 | nprim = mxalld; | |
153 | Warning("mgdraw", "nprim > mxalld, nprim: %6d pdg: %6d mreg %6d p %13.3f step %13.3f\n", | |
154 | ALLDLT.nalldl, | |
155 | fluka->PDGFromId(TRACKR.jtrack), | |
156 | mreg, | |
157 | TRACKR.ptrack, | |
158 | TRACKR.ctrack); | |
159 | } | |
160 | fluka->PrimaryIonisationStepping(nprim); | |
161 | } // primary ionisation switched on | |
162 | } // tracking resumed | |
163 | } // end of mgdraw | |
164 | } // end of extern "C" | |
165 |