PrimaryElectronIndex corrected.
[u/mrichter/AliRoot.git] / TFluka / mgdraw.cxx
CommitLineData
fa3d1cc7 1#include <Riostream.h>
8fd1d27e 2#include "TVirtualMCApplication.h"
b1b2005d 3#include "TVirtualMCStack.h"
a7bb59a2 4
fa3d1cc7 5#include "TFluka.h"
d566901f 6#include "TFlukaCodes.h"
b1b2005d 7// Fluka include
8#include "Fdimpar.h" //(DIMPAR) fluka include
9#include "Fdblprc.h" //(DBLPRC) fluka common
10#include "Ftrackr.h" //(TRACKR) fluka common
3a625972 11#include "Fopphst.h" //(OPPHST) fluka common
81f1d030 12#include "Fflkstk.h" //(FLKSTK) fluka common
d566901f 13#include "Fltclcm.h" //(LTCLCM) fluka common
14#include "Fpaprop.h" //(PAPROP) fluka common
875111c2 15#include "Falldlt.h" //(ALLDLT) fluka common
b1b2005d 16
fa3d1cc7 17#ifndef WIN32
18# define mgdraw mgdraw_
19#else
20# define mgdraw MGDRAW
21#endif
22
4aba9d66 23
24#include "TGeoManager.h" // <- delete
25
fa3d1cc7 26extern "C" {
27void mgdraw(Int_t& icode, Int_t& mreg)
28{
b1b2005d 29 TFluka* fluka = (TFluka*) gMC;
d566901f 30 if (mreg == fluka->GetDummyRegion()) return;
b1b2005d 31//
32// Make sure that stack has currrent track Id
57dc5a4a 33//
3a625972 34 Int_t trackId = -1;
b1b2005d 35 TVirtualMCStack* cppstack = fluka->GetStack();
3a625972 36
37 if (TRACKR.jtrack == -1) {
18e0cabb 38 trackId = OPPHST.louopp[OPPHST.lstopp];
39 if (trackId == 0) {
40 trackId = FLKSTK.ispark[FLKSTK.npflka][mkbmx2-1];
41 }
3a625972 42 } else {
18e0cabb 43 trackId = TRACKR.ispusr[mkbmx2-1];
44 }
45
46 Int_t verbosityLevel = fluka->GetVerbosityLevel();
47
48 if (TRACKR.jtrack < -6) {
4aba9d66 49 // from -7 to -12 = "heavy" fragment
50 // assing parent id
18e0cabb 51 // id < -6 was skipped in stuprf => if (kpart < -6) return;
52 if (verbosityLevel >= 3) {
4aba9d66 53 cout << "mgdraw: (heavy fragment) jtrack < -6 =" << TRACKR.jtrack
18e0cabb 54 << " assign parent pdg=" << fluka->PDGFromId(TRACKR.ispusr[mkbmx2 - 3]) << endl;
55 }
56 TRACKR.jtrack = TRACKR.ispusr[mkbmx2 - 3];
3a625972 57 }
58
b1b2005d 59 cppstack->SetCurrentTrack(trackId);
60//
61//
4aba9d66 62 Int_t mlttc = TRACKR.lt1trk; // LTCLCM.mlatm1;
d566901f 63 fluka->SetMreg(mreg, mlttc);
d566901f 64 fluka->SetIcode((FlukaProcessCode_t) icode);
65 fluka->SetCaller(kMGDRAW);
57dc5a4a 66
4aba9d66 67 Int_t nodeId;
875111c2 68 Int_t volId = fluka->CurrentVolID(nodeId);
4aba9d66 69 Int_t crtlttc = gGeoManager->GetCurrentNodeId()+1;
70
4aba9d66 71 // check region lattice consistency (debug Ernesto)
72 // *****************************************************
73 if( 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
57dc5a4a 86 if (!TRACKR.ispusr[mkbmx2 - 2]) {
875111c2 87 if (verbosityLevel >= 3) {
88 cout << endl << "mgdraw: energy deposition for:" << trackId
89 << " icode=" << icode
90 << " pdg=" << fluka->PDGFromId(TRACKR.jtrack)
91 << " flukaid="<< TRACKR.jtrack
92 << " mreg=" << mreg
93 << " np =" << ALLDLT.nalldl
94 << endl;
95
96 }
97 Int_t nprim = 0;
98// printf("mgdraw %5d %5d \n", ALLDLT.nalldl, ALLDLT.lalldl );
99
100 if ((nprim = ALLDLT.nalldl) > 0) {
101 //
9c1667f7 102 // Multiple steps (primary ionisation)
103 if (nprim >= mxalld) {
104 nprim = mxalld;
105 Warning("mgdraw", "nprim > mxalld, nprim: %6d pdg: %6d mreg %6d p %13.3f step %13.3f\n",
106 ALLDLT.nalldl,
107 fluka->PDGFromId(TRACKR.jtrack),
108 mreg,
109 TRACKR.ptrack,
110 TRACKR.ctrack);
111
112 }
875111c2 113 for (Int_t i = 0; i < nprim; i++) {
875111c2 114 fluka->SetCurrentPrimaryElectronIndex(i);
9c1667f7 115 (TVirtualMCApplication::Instance())->Stepping();
875111c2 116 if (i == 0) fluka->SetTrackIsNew(kFALSE);
117 }
118 fluka->SetCurrentPrimaryElectronIndex(-1);
119 } else {
120 // Single step
121 (TVirtualMCApplication::Instance())->Stepping();
122 fluka->SetTrackIsNew(kFALSE);
123 }
124
57dc5a4a 125 } else {
18e0cabb 126 //
127 // Tracking is being resumed after secondary tracking
128 //
129 if (verbosityLevel >= 3) {
130 cout << endl << "mgdraw: resuming Stepping(): " << trackId << endl;
131 }
d566901f 132
18e0cabb 133 fluka->SetTrackIsNew(kTRUE);
134 fluka->SetCaller(kMGResumedTrack);
135 (TVirtualMCApplication::Instance())->Stepping();
5d80a015 136
18e0cabb 137 // Reset flag and stored values
138 TRACKR.ispusr[mkbmx2 - 2] = 0;
139 for (Int_t i = 0; i < 9; i++) TRACKR.spausr[i] = -1.;
5d80a015 140
57dc5a4a 141
18e0cabb 142 if (verbosityLevel >= 3) {
143 cout << endl << " !!! I am in mgdraw - first Stepping() after resume: " << icode << endl;
144 cout << " Track= " << trackId << " region = " << mreg << endl;
145 }
5d80a015 146
18e0cabb 147 fluka->SetTrackIsNew(kFALSE);
148 fluka->SetCaller(kMGDRAW);
149 (TVirtualMCApplication::Instance())->Stepping();
57dc5a4a 150 }
fa3d1cc7 151} // end of mgdraw
152} // end of extern "C"
153