No double step.
[u/mrichter/AliRoot.git] / TFluka / endraw.cxx
CommitLineData
fa3d1cc7 1#include <Riostream.h>
8fd1d27e 2#include "TVirtualMCApplication.h"
df14e7fb 3#include "TGeoMaterial.h"
4#include "TGeoManager.h"
18e0cabb 5#include <TParticle.h>
df14e7fb 6#include "TFlukaCerenkov.h"
a7bb59a2 7
fa3d1cc7 8#include "TFluka.h"
d566901f 9#include "TFlukaCodes.h"
a7bb59a2 10
7dac99f1 11#include "Fdimpar.h" //(DIMPAR) fluka include
12#include "Ftrackr.h" //(TRACKR) fluka common
d566901f 13#include "Fltclcm.h" //(LTCLCM) fluka common
14#include "Fpaprop.h" //(PAPROP) fluka common
18e0cabb 15
fa3d1cc7 16#ifndef WIN32
17# define endraw endraw_
18#else
19# define endraw ENDRAW
20#endif
21extern "C" {
22void endraw(Int_t& icode, Int_t& mreg, Double_t& rull, Double_t& xsco, Double_t& ysco, Double_t& zsco)
23{
4b81a588 24 TFluka* fluka = (TFluka*) gMC;
d566901f 25 // nothing to do if particle in dummy region
26 if (mreg == fluka->GetDummyRegion()) return;
2bc4c610 27 Int_t verbosityLevel = fluka->GetVerbosityLevel();
2bc78182 28 Bool_t debug = (verbosityLevel >= 3)? kTRUE : kFALSE;
4aba9d66 29 Int_t mlttc = (icode==kKASKADinelarecoil) ? TRACKR.lt2trk : TRACKR.lt1trk; //LTCLCM.mlatm1;
d566901f 30 fluka->SetCaller(kENDRAW);
4b81a588 31 fluka->SetRull(rull);
32 fluka->SetXsco(xsco);
33 fluka->SetYsco(ysco);
34 fluka->SetZsco(zsco);
d566901f 35 fluka->SetMreg(mreg, mlttc);
4aba9d66 36
022d4275 37 Float_t edep = rull;
69e406a7 38 Int_t ipt = fluka->PDGFromId(TRACKR.jtrack);
39 if (ipt == -1) {
40 if (debug) printf("Unknown particle %5d %5d \n", TRACKR.jtrack, icode);
41 return;
42 }
43
022d4275 44
df14e7fb 45 if (TRACKR.jtrack == -1) {
18e0cabb 46 // Handle quantum efficiency the G3 way
2bc4c610 47 if (debug) printf("endraw: Cerenkov photon depositing energy: %d %e\n", mreg, rull);
df14e7fb 48 TGeoMaterial* material = (gGeoManager->GetCurrentVolume())->GetMaterial();
df14e7fb 49 TFlukaCerenkov* cerenkov = dynamic_cast<TFlukaCerenkov*> (material->GetCerenkovProperties());
50 if (cerenkov) {
18e0cabb 51 Double_t eff = (cerenkov->GetQuantumEfficiency(rull));
52 if (gRandom->Rndm() > eff) {
53 edep = 0.;
54 }
df14e7fb 55 }
56 }
2bc78182 57
18e0cabb 58 TVirtualMCStack* cppstack = fluka->GetStack();
18e0cabb 59
60 if (debug) {
4aba9d66 61 cout << "ENDRAW For icode=" << icode
18e0cabb 62 << " track=" << TRACKR.ispusr[mkbmx2-1] << " pdg=" << fluka->PDGFromId(TRACKR.jtrack)
4aba9d66 63 << " edep="<< edep << " mreg=" << mreg << endl;
18e0cabb 64 }
65
4aba9d66 66 // check region lattice consistency (debug Ernesto)
67 // *****************************************************
68 Int_t nodeId;
69 Int_t volId = fluka->CurrentVolID(nodeId);
70 Int_t crtlttc = gGeoManager->GetCurrentNodeId()+1;
4b1c3976 71 if(debug && mreg != volId && !gGeoManager->IsOutside() ) {
4aba9d66 72 cout << " endraw: track=" << TRACKR.ispusr[mkbmx2-1] << " pdg=" << fluka->PDGFromId(TRACKR.jtrack)
73 << " icode=" << icode << " gNstep=" << fluka->GetNstep() << endl
74 << " fluka mreg=" << mreg << " mlttc=" << mlttc << endl
75 << " TGeo volId=" << volId << " crtlttc=" << crtlttc << endl
76 << " common TRACKR lt1trk=" << TRACKR.lt1trk << " lt2trk=" << TRACKR.lt2trk << endl
77 << " common LTCLCM newlat=" << LTCLCM.newlat << " mlatld=" << LTCLCM.mlatld << endl
78 << " mlatm1=" << LTCLCM.mlatm1 << " mltsen=" << LTCLCM.mltsen << endl
79 << " mltsm1=" << LTCLCM.mltsm1 << " mlattc=" << LTCLCM.mlattc << endl;
80 if( mlttc == crtlttc ) cout << " *************************************************************" << endl;
81 }
82 // *****************************************************
d566901f 83 if (icode != kEMFSCOstopping1 && icode != kEMFSCOstopping2) {
84 fluka->SetIcode((FlukaProcessCode_t)icode);
022d4275 85 fluka->SetRull(edep);
18e0cabb 86 if (icode == kKASKADelarecoil && TRACKR.ispusr[mkbmx2-5]) {
69e406a7 87 // Elastic recoil and in stuprf npprmr > 0,
88 // the secondary being loaded is actually still the interacting particle
89 cppstack->SetCurrentTrack( TRACKR.ispusr[mkbmx2-4] );
90 } else {
4aba9d66 91 cppstack->SetCurrentTrack(TRACKR.ispusr[mkbmx2-1] );
69e406a7 92 }
93
95bb1461 94 (TVirtualMCApplication::Instance())->Stepping();
69e406a7 95
022d4275 96 } else {
2bc78182 97 //
d566901f 98 // For icode 21,22 the particle has fallen below thresshold.
2bc78182 99 // This has to be signalled to the StepManager()
100 //
18e0cabb 101 cppstack->SetCurrentTrack( TRACKR.ispusr[mkbmx2-1] );
a5b73d69 102 fluka->SetRull(edep);
d566901f 103 fluka->SetIcode((FlukaProcessCode_t) icode);
022d4275 104 (TVirtualMCApplication::Instance())->Stepping();
105 fluka->SetTrackIsNew(kFALSE);
0375c7c5 106// fluka->SetIcode((FlukaProcessCode_t)icode);
107// fluka->SetRull(0.);
108// (TVirtualMCApplication::Instance())->Stepping();
022d4275 109 }
fa3d1cc7 110} // end of endraw
111} // end of extern "C"
112