]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TFluka/stupre.cxx
Pemf generation option is obsolete. Removed.
[u/mrichter/AliRoot.git] / TFluka / stupre.cxx
index 9183498d3b2ba349eeceb13ac518a429084ec1a8..b1f54559fecd08653a7de3fafa07a80a0b8d0158 100644 (file)
@@ -18,6 +18,7 @@
 
 //Virtual MC
 #include "TFluka.h"
+#include "TFlukaCodes.h"
 #include "TVirtualMCStack.h"
 #include "TVirtualMCApplication.h"
 #include "TParticle.h"
@@ -31,13 +32,15 @@ void stupre()
 //*  SeT User PRoperties for Emf particles                               *
 //*                                                                      *
 //*----------------------------------------------------------------------*
-  static Double_t emassmev = PAPROP.am[9] * 1000.;
+// Get the pointer to the VMC
+  TFluka* fluka =  (TFluka*) gMC;
+  static Double_t emassmev = PAPROP.am[kFLUKAelectron + 6] * 1000.;
     
   Int_t lbhabh = 0;
   if (EVTFLG.ldltry == 1) {
     if (EMFSTK.ichemf[EMFSTK.npemf-1] * EMFSTK.ichemf[EMFSTK.npemf-2] < 0) lbhabh = 1;
   }
-  
+
 // mkbmx1         = dimension for kwb real spare array in fluka stack in DIMPAR
 // mkbmx2         = dimension for kwb int. spare array in fluka stack in DIMPAR
 // EMFSTK.espark  = spare real variables available for 
@@ -50,17 +53,16 @@ void stupre()
   Int_t npnw, ispr;
   for (npnw = EMFSTK.npstrt-1; npnw <= EMFSTK.npemf-1; npnw++) {
 
-    for (ispr = 0; ispr <= mkbmx1-1; ispr++) 
-      EMFSTK.espark[npnw][ispr] = TRACKR.spausr[ispr];
-    
-    for (ispr = 0; ispr <= mkbmx2-1; ispr++) 
-      EMFSTK.iespak[npnw][ispr] = TRACKR.ispusr[ispr];
-    
-    EMFSTK.louemf[npnw] = TRACKR.llouse;
+      for (ispr = 0; ispr <= mkbmx1-1; ispr++) 
+          EMFSTK.espark[npnw][ispr] = TRACKR.spausr[ispr];
+
+      for (ispr = 0; ispr <= mkbmx2-1; ispr++) 
+          EMFSTK.iespak[npnw][ispr] = TRACKR.ispusr[ispr];
+
+      EMFSTK.louemf[npnw] = TRACKR.llouse;
   }
 
-// Get the pointer to the VMC
-  TFluka* fluka =  (TFluka*) gMC;
+
   Int_t verbosityLevel = fluka->GetVerbosityLevel();
   Bool_t debug = (verbosityLevel>=3)?kTRUE:kFALSE;
   fluka->SetTrackIsNew(kTRUE);
@@ -75,40 +77,43 @@ void stupre()
   
   for (kp = EMFSTK.npstrt - 1; kp <= EMFSTK.npemf - 1; kp++) {
     
-// Ckeck transport cut first
-    Int_t ireg = EMFSTK.iremf[kp];
-    Double_t cut = (TMath::Abs(EMFSTK.ichemf[kp]) == 1) ? EFMRGN.elethr[ireg-1] :  EFMRGN.phothr[ireg-1];
+// Check transport cut first
+    Int_t    ireg   = EMFSTK.iremf[kp];
+    Double_t cut    = (TMath::Abs(EMFSTK.ichemf[kp]) == 1) ? EMFRGN.elethr[ireg-1] :  EMFRGN.phothr[ireg-1];
     Double_t e      = EMFSTK.etemf[kp];
+
     if ((e < cut) 
-       && ( 
-           (EMFSTK.ichemf[kp] ==  0) ||
-           (EMFSTK.ichemf[kp] == -1) ||
-           (EMFSTK.ichemf[kp] ==  1 &&  EFMRGN.phothr[ireg-1] > emassmev)
-           )
-       )
+        && (
+            (EMFSTK.ichemf[kp] ==  0) ||
+            (EMFSTK.ichemf[kp] == -1) ||
+            (EMFSTK.ichemf[kp] ==  1 &&  EMFRGN.phothr[ireg-1] > emassmev)
+            )
+        )
     {
-       EMFSTK.iespak[kp][mkbmx2-1] = -1;
-       EMFSTK.iespak[kp][mkbmx2-2] =  0;
-       continue;
+        EMFSTK.iespak[kp][mkbmx2-1] = TRACKR.ispusr[mkbmx2-1];
+        EMFSTK.iespak[kp][mkbmx2-2] =  0;
+        EMFSTK.iespak[kp][mkbmx2-3] = TRACKR.jtrack;
+        EMFSTK.irlatt[kp] = TRACKR.lt1trk;
+        continue;
     }
 
-//* save the parent track number and reset it at each loop
+// Save the parent track number and reset it at each loop
     Int_t done = 0;
     Int_t parent =  TRACKR.ispusr[mkbmx2-1];
     Int_t flukaid = 0;
-    
+
 // Identify particle type
-    if      (EMFSTK.ichemf[kp] == -1)  flukaid = 3;
-    else if (EMFSTK.ichemf[kp] ==  0)  flukaid = 7;
-    else if (EMFSTK.ichemf[kp] ==  1)  flukaid = 4;
-    
+    if      (EMFSTK.ichemf[kp] == -1)  flukaid = kFLUKAelectron;
+    else if (EMFSTK.ichemf[kp] ==  0)  flukaid = kFLUKAphoton;
+    else if (EMFSTK.ichemf[kp] ==  1)  flukaid = kFLUKApositron;
 
-    e      *= emvgev;
+
+    e *= emvgev;
     Int_t    pdg    = fluka->PDGFromId(flukaid);
     Double_t p      = sqrt(e * e - PAPROP.am[flukaid+6] * PAPROP.am[flukaid+6]);
     Double_t px     = p * EMFSTK.uemf[kp];
-    Double_t pz     = p * EMFSTK.vemf[kp];
-    Double_t py     = p * EMFSTK.wemf[kp];
+    Double_t py     = p * EMFSTK.vemf[kp];
+    Double_t pz     = p * EMFSTK.wemf[kp];
     Double_t tof    = EMFSTK.agemf[kp];
     Double_t polx   = EMFSTK.upol[kp];
     Double_t poly   = EMFSTK.vpol[kp];
@@ -126,73 +131,73 @@ void stupre()
 //* all secondaries are true
     if ((EVTFLG.lpairp == 1) || (EVTFLG.lphoel == 1) ||
         (EVTFLG.lannfl == 1) || (EVTFLG.lannrs == 1)) {
-       
-       if (EVTFLG.lpairp == 1) mech = kPPair;
-       else if (EVTFLG.lphoel == 1) mech = kPPhotoelectric;
-       else mech = kPAnnihilation;
+
+        if (EVTFLG.lpairp == 1) mech = kPPair;
+        else if (EVTFLG.lphoel == 1) mech = kPPhotoelectric;
+        else mech = kPAnnihilation;
         cppstack->PushTrack(done, parent, pdg,
-                          px, py, pz, e, vx, vy, vz, tof,
-                          polx, poly, polz, mech, ntr, weight, is);
-       if (debug) cout << endl << " !!! stupre (PAIR, ..) : ntr=" << ntr << "pdg " << pdg << " parent=" << parent << "energy " << e-PAPROP.am[flukaid+6] << endl;
+                           px, py, pz, e, vx, vy, vz, tof,
+                           polx, poly, polz, mech, ntr, weight, is);
+        if (debug) cout << endl << " !!! stupre (PAIR, ..) : ntr=" << ntr << " pdg " << pdg << " parent=" << parent << " energy " << e-PAPROP.am[flukaid+6] << endl;
 
-       EMFSTK.iespak[kp][mkbmx2-1] = ntr;
-       EMFSTK.iespak[kp][mkbmx2-2] = 0;
+        EMFSTK.iespak[kp][mkbmx2-1] = ntr;
+        EMFSTK.iespak[kp][mkbmx2-2] = 0;
     } // end of lpairp, lphoel, lannfl, lannrs
-    
+
 //* Compton: secondary is true only if charged (e+, e-)
     else if ((EVTFLG.lcmptn == 1)) {
 
-       if (EMFSTK.ichemf[kp] != 0) {
-           mech = kPCompton;
-           cppstack->PushTrack(done, parent, pdg,
-                              px, py, pz, e, vx, vy, vz, tof,
-                              polx, poly, polz, mech, ntr, weight, is);
-           if (debug) cout << endl << " !!! stupre (COMPTON) : ntr=" << ntr << "pdg " << pdg << " parent=" << parent << endl;
-           EMFSTK.iespak[kp][mkbmx2-1] = ntr;
-           EMFSTK.iespak[kp][mkbmx2-2] = 0;
-       }
+        if (EMFSTK.ichemf[kp] != 0) {
+            mech = kPCompton;
+            cppstack->PushTrack(done, parent, pdg,
+                               px, py, pz, e, vx, vy, vz, tof,
+                               polx, poly, polz, mech, ntr, weight, is);
+            if (debug) cout << endl << " !!! stupre (COMPTON) : ntr=" << ntr << " pdg " << pdg << " parent=" << parent << endl;
+            EMFSTK.iespak[kp][mkbmx2-1] = ntr;
+            EMFSTK.iespak[kp][mkbmx2-2] = 0;
+        }
     } // end of lcmptn
-    
+
 //* Bremsstrahlung: true secondary only if charge = 0 (photon)
     else if ((EVTFLG.lbrmsp == 1)) {
-       if (EMFSTK.ichemf[kp] == 0) {
-           mech = kPBrem;
-           cppstack->PushTrack(done, parent, pdg,
-                              px, py, pz, e, vx, vy, vz, tof,
-                              polx, poly, polz, mech, ntr, weight, is);
-           if (debug) cout << endl << " !!! stupre (BREMS) : ntr=" << ntr << "pdg " << pdg << " parent=" << parent << endl;
-           EMFSTK.iespak[kp][mkbmx2-1] = ntr;
-           EMFSTK.iespak[kp][mkbmx2-2] = 0;
-       }
+        if (EMFSTK.ichemf[kp] == 0) {
+            mech = kPBrem;
+            cppstack->PushTrack(done, parent, pdg,
+                               px, py, pz, e, vx, vy, vz, tof,
+                               polx, poly, polz, mech, ntr, weight, is);
+            if (debug) cout << endl << " !!! stupre (BREMS) : ntr=" << ntr << " pdg " << pdg << " parent=" << parent << endl;
+            EMFSTK.iespak[kp][mkbmx2-1] = ntr;
+            EMFSTK.iespak[kp][mkbmx2-2] = 0;
+        }
     } // end of lbrmsp
-    
+
 //* Delta ray: If Bhabha, true secondary only if negative (electron)
     else if ((EVTFLG.ldltry == 1)) {
-       if (lbhabh == 1) {
-           if (EMFSTK.ichemf[kp] == -1) {
-               mech = kPDeltaRay;
-               cppstack->PushTrack(done, parent, pdg,
-                                  px, py, pz, e, vx, vy, vz, tof,
-                                  polx, poly, polz, mech, ntr, weight, is);
-               EMFSTK.iespak[kp][mkbmx2-1] = ntr;
-               EMFSTK.iespak[kp][mkbmx2-2] = 0;
-          if (debug) cout << endl << " !!! stupre (BHABA) : ntr=" << ntr << "pdg " << pdg << " parent=" << parent << endl;
-           } // end of Bhabha
-       } // lbhabh == 1
-       
+        if (lbhabh == 1) {
+            if (EMFSTK.ichemf[kp] == -1) {
+                mech = kPDeltaRay;
+                cppstack->PushTrack(done, parent, pdg,
+                                   px, py, pz, e, vx, vy, vz, tof,
+                                   polx, poly, polz, mech, ntr, weight, is);
+                EMFSTK.iespak[kp][mkbmx2-1] = ntr;
+                EMFSTK.iespak[kp][mkbmx2-2] = 0;
+           if (debug) cout << endl << " !!! stupre (BHABA) : ntr=" << ntr << " pdg " << pdg << " parent=" << parent << endl;
+            } // end of Bhabha
+        } // lbhabh == 1
+
 //* Delta ray: Otherwise Moller: true secondary is the electron with
 //*            lower energy, which has been put higher in the stack
-       else if (kp == EMFSTK.npemf-1) {
-           mech = kPDeltaRay;
-           cppstack->PushTrack(done, parent, pdg,
-                              px, py, pz, e, vx, vy, vz, tof,
-                              polx, poly, polz, mech, ntr, weight, is);
-           if (debug) cout << endl << " !!! stupre (Moller) : ntr=" << ntr << "pdg " << pdg << " parent=" << parent << endl;
-           EMFSTK.iespak[kp][mkbmx2-1] = ntr;
-           EMFSTK.iespak[kp][mkbmx2-2] = 0;
-       } // end of Delta ray
+        else if (kp == EMFSTK.npemf-1) {
+            mech = kPDeltaRay;
+            cppstack->PushTrack(done, parent, pdg,
+                               px, py, pz, e, vx, vy, vz, tof,
+                               polx, poly, polz, mech, ntr, weight, is);
+            if (debug) cout << endl << " !!! stupre (Moller) : ntr=" << ntr << " pdg " << pdg << " parent=" << parent << endl;
+            EMFSTK.iespak[kp][mkbmx2-1] = ntr;
+            EMFSTK.iespak[kp][mkbmx2-2] = 0;
+        } // end of Delta ray
     } // end of ldltry
-    
+
   } // end of loop
 } // end of stupre
 } // end of extern "C"