Updated buspatch and DDL numbers for station 345 and started buspatch at 1
[u/mrichter/AliRoot.git] / TFluka / mgdraw.cxx
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
49     if (TRACKR.jtrack < -6) {
50        // from -7 to -12 = "heavy" fragment
51        // assing parent id
52        // id < -6 was skipped in stuprf =>   if (kpart < -6) return;
53        if (verbosityLevel >= 3) {
54           cout << "mgdraw: (heavy fragment) jtrack < -6 =" << TRACKR.jtrack
55                << " assign parent pdg=" << fluka->PDGFromId(TRACKR.ispusr[mkbmx2 - 3]) << endl;
56        }
57        TRACKR.jtrack = TRACKR.ispusr[mkbmx2 - 3];
58     }
59     
60     cppstack->SetCurrentTrack(trackId);
61 //
62 //    
63     Int_t mlttc = TRACKR.lt1trk; // LTCLCM.mlatm1;
64     fluka->SetMreg(mreg, mlttc);
65     fluka->SetIcode((FlukaProcessCode_t) icode);
66     fluka->SetCaller(kMGDRAW);
67
68     Int_t nodeId;
69     Int_t volId   = fluka->CurrentVolID(nodeId);
70     Int_t crtlttc = gGeoManager->GetCurrentNodeId()+1;
71
72     // check region lattice consistency (debug Ernesto)
73     // *****************************************************
74     if(verbosityLevel>=3 && mreg != volId  && !gGeoManager->IsOutside() ) {
75        cout << "  mgdraw:   track=" << trackId << " pdg=" << fluka->PDGFromId(TRACKR.jtrack)
76             << " icode=" << icode << " gNstep=" << fluka->GetNstep() << endl
77             << "               fluka   mreg=" << mreg << " mlttc=" << mlttc << endl
78             << "               TGeo   volId=" << volId << " crtlttc=" << crtlttc << endl
79             << "     common TRACKR   lt1trk=" << TRACKR.lt1trk << " lt2trk=" << TRACKR.lt2trk << endl
80             << "     common LTCLCM   newlat=" << LTCLCM.newlat << " mlatld=" <<  LTCLCM.mlatld << endl
81             << "                     mlatm1=" << LTCLCM.mlatm1 << " mltsen=" <<  LTCLCM.mltsen << endl
82             << "                     mltsm1=" << LTCLCM.mltsm1 << " mlattc=" << LTCLCM.mlattc << endl;
83         if( mlttc == crtlttc ) cout << "   *************************************************************" << endl;
84     }
85     // *****************************************************
86
87     
88     Int_t med   = FLKMAT.medium[mreg - 1];  // Medium
89     Int_t msd   = DPDXCM.msdpdx[med  - 1];  // Iionisation model
90
91     if (!TRACKR.ispusr[mkbmx2 - 2]) {
92         if (verbosityLevel >= 3) {
93             cout << endl << "mgdraw: energy deposition for:" << trackId
94                  << " icode=" << icode
95                  << " pdg=" << fluka->PDGFromId(TRACKR.jtrack)
96                  << " flukaid="<< TRACKR.jtrack
97                  << " mreg=" << mreg
98                  << " np  =" << ALLDLT.nalldl
99                  << endl;
100                  
101         }
102
103         if (msd > 0) {
104 //
105 // Primary ionsiations
106             Int_t nprim = ALLDLT.nalldl;
107 // Protection against nprim > mxalld
108             if (nprim >= mxalld) {
109                 nprim = mxalld;
110                 Warning("mgdraw", "nprim > mxalld, nprim: %6d  pdg: %6d mreg %6d p %13.3f step %13.3f\n", 
111                        ALLDLT.nalldl, 
112                        fluka->PDGFromId(TRACKR.jtrack), 
113                        mreg, 
114                        TRACKR.ptrack, 
115                        TRACKR.ctrack);
116                 
117             }
118 // Multiple steps for nprim > 0
119             if (nprim > 0) {
120                 for (Int_t i = 0; i < nprim; i++) {
121                     fluka->SetCurrentPrimaryElectronIndex(i);
122                     (TVirtualMCApplication::Instance())->Stepping();
123                     if (i == 0) fluka->SetTrackIsNew(kFALSE);
124                 }
125             } else {
126                 // No primary electron ionisation
127                 // Call Stepping anyway but flag nprim = 0 as index = -2
128                 fluka->SetCurrentPrimaryElectronIndex(-2);
129                 (TVirtualMCApplication::Instance())->Stepping();
130             }
131             // Reset the index
132             fluka->SetCurrentPrimaryElectronIndex(-1);
133         } else {
134             // Single step
135             (TVirtualMCApplication::Instance())->Stepping();
136             fluka->SetTrackIsNew(kFALSE);
137         }
138         
139     } else {
140         //
141         // Tracking is being resumed after secondary tracking
142         //
143         if (verbosityLevel >= 3) {
144             cout << endl << "mgdraw: resuming Stepping(): " << trackId << endl;
145         }
146
147         fluka->SetTrackIsNew(kTRUE);
148         fluka->SetCaller(kMGResumedTrack);
149         (TVirtualMCApplication::Instance())->Stepping();
150
151         // Reset flag and stored values
152         TRACKR.ispusr[mkbmx2 - 2] = 0;
153         for (Int_t i = 0; i < 9; i++) TRACKR.spausr[i] = -1.;
154
155
156         if (verbosityLevel >= 3) {
157             cout << endl << " !!! I am in mgdraw - first Stepping() after resume: " << icode << endl;
158             cout << " Track= " << trackId << " region = " << mreg << endl;
159         }
160
161         fluka->SetTrackIsNew(kFALSE);
162         fluka->SetCaller(kMGDRAW);
163         if (msd > 0) fluka->SetCurrentPrimaryElectronIndex(-2);
164         (TVirtualMCApplication::Instance())->Stepping();
165         if (msd > 0) fluka->SetCurrentPrimaryElectronIndex(-1);
166     }
167 } // end of mgdraw
168 } // end of extern "C"
169