Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/standard/include/G4GoudsmitSaundersonMscModel.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /processes/electromagnetic/standard/include/G4GoudsmitSaundersonMscModel.hh (Version 11.3.0) and /processes/electromagnetic/standard/include/G4GoudsmitSaundersonMscModel.hh (Version 10.6.p3)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                         3 // * License and Disclaimer                                           *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th      5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided      6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License      7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/      8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.           9 // * include a list of copyright holders.                             *
 10 // *                                               10 // *                                                                  *
 11 // * Neither the authors of this software syst     11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file      15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitatio     16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                               17 // *                                                                  *
 18 // * This  code  implementation is the result      18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboratio     19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distri     20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  ag     21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publicati     22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Sof     23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // *******************************************     24 // ********************************************************************
 25 //                                                 25 //
 26 //                                                 26 //
 27 // -------------------------------------------     27 // ----------------------------------------------------------------------------
 28 //                                                 28 //
 29 // GEANT4 Class header file                        29 // GEANT4 Class header file
 30 //                                                 30 //
 31 // File name:     G4GoudsmitSaundersonMscModel     31 // File name:     G4GoudsmitSaundersonMscModel
 32 //                                                 32 //
 33 // Author:        Mihaly Novak / (Omrane Kadri     33 // Author:        Mihaly Novak / (Omrane Kadri)
 34 //                                                 34 //
 35 // Creation date: 20.02.2009                       35 // Creation date: 20.02.2009
 36 //                                                 36 //
 37 // Modifications:                                  37 // Modifications:
 38 // 04.03.2009 V.Ivanchenko cleanup and format      38 // 04.03.2009 V.Ivanchenko cleanup and format according to Geant4 EM style
 39 // 12.05.2010 O.Kadri: adding Qn1 and Qn12 as      39 // 12.05.2010 O.Kadri: adding Qn1 and Qn12 as private doubles
 40 // 18.05.2015 M. Novak provide PLERIMINARYY ve     40 // 18.05.2015 M. Novak provide PLERIMINARYY version of updated class.
 41 //            All algorithms of the class were     41 //            All algorithms of the class were revised and updated, new methods added.
 42 //            A new version of Kawrakow-Bielaj     42 //            A new version of Kawrakow-Bielajew Goudsmit-Saunderson MSC model
 43 //            based on the screened Rutherford     43 //            based on the screened Rutherford DCS for elastic scattering of
 44 //            electrons/positrons has been int     44 //            electrons/positrons has been introduced[1,2]. The corresponding MSC
 45 //            angular distributions over a 2D      45 //            angular distributions over a 2D parameter grid have been recomputed
 46 //            and the CDFs are now stored in a     46 //            and the CDFs are now stored in a variable transformed (smooth) form[2,3]
 47 //            together with the corresponding      47 //            together with the corresponding rational interpolation parameters.
 48 //            These angular distributions are      48 //            These angular distributions are handled by the new
 49 //            G4GoudsmitSaundersonTable class      49 //            G4GoudsmitSaundersonTable class that is responsible to sample if
 50 //            it was no, single, few or multip     50 //            it was no, single, few or multiple scattering case and delivers the
 51 //            angular deflection (i.e. cos(the     51 //            angular deflection (i.e. cos(theta) and sin(theta)).
 52 //            Two screening options are provid     52 //            Two screening options are provided:
 53 //             - if fgIsUsePWATotalXsecData=TR     53 //             - if fgIsUsePWATotalXsecData=TRUE i.e. SetOptionPWAScreening(TRUE)
 54 //               was called before initialisat     54 //               was called before initialisation: screening parameter value A is
 55 //               determined such that the firs     55 //               determined such that the first transport coefficient G1(A)
 56 //               computed according to the scr     56 //               computed according to the screened Rutherford DCS for elastic
 57 //               scattering will reproduce the     57 //               scattering will reproduce the one computed from the PWA elastic
 58 //               and first transport mean free     58 //               and first transport mean free paths[4].
 59 //             - if fgIsUsePWATotalXsecData=FA     59 //             - if fgIsUsePWATotalXsecData=FALSE i.e. default value or
 60 //               SetOptionPWAScreening(FALSE)      60 //               SetOptionPWAScreening(FALSE) was called before initialisation:
 61 //               screening parameter value A i     61 //               screening parameter value A is computed according to Moliere's
 62 //               formula (by using material de     62 //               formula (by using material dependent parameters \chi_cc2 and b_c
 63 //               precomputed for each material     63 //               precomputed for each material used at initialization in
 64 //               G4GoudsmitSaundersonTable) [3     64 //               G4GoudsmitSaundersonTable) [3]
 65 //            Elastic and first trasport mean      65 //            Elastic and first trasport mean free paths are used consistently.
 66 //            The new version is self-consiste     66 //            The new version is self-consistent, several times faster, more
 67 //            robust and accurate compared to      67 //            robust and accurate compared to the earlier version.
 68 //            Spin effects as well as a more a     68 //            Spin effects as well as a more accurate energy loss correction and
 69 //            computations of Lewis moments wi     69 //            computations of Lewis moments will be implemented later on.
 70 // 02.09.2015 M. Novak: first version of new s     70 // 02.09.2015 M. Novak: first version of new step limit is provided.
 71 //            fUseSafetyPlus corresponds to Ur     71 //            fUseSafetyPlus corresponds to Urban fUseSafety (default)
 72 //            fUseDistanceToBoundary correspon     72 //            fUseDistanceToBoundary corresponds to Urban fUseDistanceToBoundary
 73 //            fUseSafety  corresponds to EGSnr     73 //            fUseSafety  corresponds to EGSnrc error-free stepping algorithm
 74 //            Range factor can be significantl     74 //            Range factor can be significantly higher at each case than in Urban.
 75 // 23.08.2017 M. Novak: added corrections to a     75 // 23.08.2017 M. Novak: added corrections to account spin effects (Mott-correction).
 76 //            It can be activated by setting t     76 //            It can be activated by setting the fIsMottCorrection flag to be true
 77 //            before initialization using the      77 //            before initialization using the SetOptionMottCorrection() public method.
 78 //            The fMottCorrection member is re     78 //            The fMottCorrection member is responsible to handle pre-computed Mott
 79 //            correction (rejection) functions     79 //            correction (rejection) functions obtained by numerically computing
 80 //            Goudsmit-Saunderson agnular dist     80 //            Goudsmit-Saunderson agnular distributions based on a DCS accounting spin
 81 //            effects and screening correction     81 //            effects and screening corrections. The DCS used to compute the accurate
 82 //            GS angular distributions is: DCS     82 //            GS angular distributions is: DCS_{cor} = DCS_{SR}x[ DCS_{R}/DCS_{Mott}] where :
 83 //               # DCS_{SR} is the relativisti     83 //               # DCS_{SR} is the relativistic Screened-Rutherford DCS (first Born approximate
 84 //                 solution of the Klein-Gordo     84 //                 solution of the Klein-Gordon i.e. relativistic Schrodinger equation =>
 85 //                 scattering of spinless e- o     85 //                 scattering of spinless e- on exponentially screened Coulomb potential)
 86 //                 note: the default (without      86 //                 note: the default (without using Mott-correction) GS angular distributions
 87 //                 are based on this DCS_{SR}      87 //                 are based on this DCS_{SR} with Moliere's screening parameter!
 88 //               # DCS_{R} is the Rutherford D     88 //               # DCS_{R} is the Rutherford DCS which is the same as above but without
 89 //                 screening                       89 //                 screening
 90 //               # DCS_{Mott} is the Mott DCS      90 //               # DCS_{Mott} is the Mott DCS i.e. solution of the Dirac equation with a bare
 91 //                 Coulomb potential i.e. scat     91 //                 Coulomb potential i.e. scattering of particles with spin (e- or e+) on a
 92 //                 point-like unscreened Coulo     92 //                 point-like unscreened Coulomb potential
 93 //               # moreover, the screening par     93 //               # moreover, the screening parameter of the DCS_{cor} was determined such that
 94 //                 the DCS_{cor} with this cor     94 //                 the DCS_{cor} with this corrected screening parameter reproduce the first
 95 //                 transport cross sections ob     95 //                 transport cross sections obtained from the corresponding most accurate DCS
 96 //                 (i.e. from elsepa [4])          96 //                 (i.e. from elsepa [4])
 97 //            Unlike the default GS, the Mott-     97 //            Unlike the default GS, the Mott-corrected angular distributions are particle type
 98 //            (different for e- and e+ <= the      98 //            (different for e- and e+ <= the DCS_{Mott} and the screening correction) and target
 99 //            (Z and material) dependent.          99 //            (Z and material) dependent.
100 // 02.02.2018 M. Novak: implemented CrossSecti    100 // 02.02.2018 M. Novak: implemented CrossSectionPerVolume interface method (used only for testing)
101 //                                                101 //
102 // Class description:                             102 // Class description:
103 //   Kawrakow-Bielajew Goudsmit-Saunderson MSC    103 //   Kawrakow-Bielajew Goudsmit-Saunderson MSC model based on the screened Rutherford DCS
104 //   for elastic scattering of e-/e+. Option,     104 //   for elastic scattering of e-/e+. Option, to include (Mott) correction (see above), is
105 //   also available now (SetOptionMottCorrecti    105 //   also available now (SetOptionMottCorrection(true)). An EGSnrc like error-free stepping
106 //   algorithm (UseSafety) is available beyond    106 //   algorithm (UseSafety) is available beyond the usual Geant4 step limitation algorithms
107 //   and true to geomerty and geometry to true    107 //   and true to geomerty and geometry to true step length computations that were adopted
108 //   from the Urban model[5]. The most accurat    108 //   from the Urban model[5]. The most accurate setting: error-free stepping (UseSafety)
109 //   with Mott-correction (SetOptionMottCorrec    109 //   with Mott-correction (SetOptionMottCorrection(true)).
110 //                                                110 //
111 // References:                                    111 // References:
112 //   [1] A.F.Bielajew, NIMB 111 (1996) 195-208    112 //   [1] A.F.Bielajew, NIMB 111 (1996) 195-208
113 //   [2] I.Kawrakow, A.F.Bielajew, NIMB 134(19    113 //   [2] I.Kawrakow, A.F.Bielajew, NIMB 134(1998) 325-336
114 //   [3] I.Kawrakow, E.Mainegra-Hing, D.W.O.Ro    114 //   [3] I.Kawrakow, E.Mainegra-Hing, D.W.O.Rogers, F.Tessier,B.R.B.Walters, NRCC
115 //       Report PIRS-701 (2013)                   115 //       Report PIRS-701 (2013)
116 //   [4] F.Salvat, A.Jablonski, C.J. Powell, C    116 //   [4] F.Salvat, A.Jablonski, C.J. Powell, CPC 165(2005) 157-190
117 //   [5] L.Urban, Preprint CERN-OPEN-2006-077     117 //   [5] L.Urban, Preprint CERN-OPEN-2006-077 (2006)
118 //                                                118 //
119 // -------------------------------------------    119 // -----------------------------------------------------------------------------
120                                                   120 
121 #ifndef G4GoudsmitSaundersonMscModel_h            121 #ifndef G4GoudsmitSaundersonMscModel_h
122 #define G4GoudsmitSaundersonMscModel_h 1          122 #define G4GoudsmitSaundersonMscModel_h 1
123                                                   123 
124 #include <CLHEP/Units/SystemOfUnits.h>            124 #include <CLHEP/Units/SystemOfUnits.h>
125                                                   125 
126 #include "G4VMscModel.hh"                         126 #include "G4VMscModel.hh"
127 #include "G4PhysicsTable.hh"                      127 #include "G4PhysicsTable.hh"
128 #include "G4MaterialCutsCouple.hh"                128 #include "G4MaterialCutsCouple.hh"
129 #include "globals.hh"                             129 #include "globals.hh"
130                                                   130 
131                                                   131 
132 class G4DataVector;                               132 class G4DataVector;
133 class G4ParticleChangeForMSC;                     133 class G4ParticleChangeForMSC;
134 class G4LossTableManager;                         134 class G4LossTableManager;
135 class G4GoudsmitSaundersonTable;                  135 class G4GoudsmitSaundersonTable;
136 class G4GSPWACorrections;                         136 class G4GSPWACorrections;
137                                                   137 
138 class G4GoudsmitSaundersonMscModel : public G4    138 class G4GoudsmitSaundersonMscModel : public G4VMscModel
139 {                                                 139 {
140 public:                                           140 public:
141                                                   141 
142   G4GoudsmitSaundersonMscModel(const G4String&    142   G4GoudsmitSaundersonMscModel(const G4String& nam = "GoudsmitSaunderson");
143                                                   143 
144   ~G4GoudsmitSaundersonMscModel() override;    << 144   virtual ~G4GoudsmitSaundersonMscModel();
145                                                   145 
146   void Initialise(const G4ParticleDefinition*, << 146   virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&);
147                                                   147 
148   void InitialiseLocal(const G4ParticleDefinit << 148   virtual void InitialiseLocal(const G4ParticleDefinition* p, G4VEmModel* masterModel);
149                                                   149 
150   G4ThreeVector& SampleScattering(const G4Thre << 
151                                                   150 
152   G4double ComputeTruePathLengthLimit(const G4 << 151   virtual G4ThreeVector& SampleScattering(const G4ThreeVector&, G4double safety);
153                                                   152 
154   G4double ComputeGeomPathLength(G4double true << 153   virtual G4double ComputeTruePathLengthLimit(const G4Track& track, G4double& currentMinimalStep);
155                                                   154 
156   G4double ComputeTrueStepLength(G4double geom << 155   virtual G4double ComputeGeomPathLength(G4double truePathLength);
                                                   >> 156 
                                                   >> 157   virtual G4double ComputeTrueStepLength(G4double geomStepLength);
157                                                   158 
158   // method to compute first transport cross s    159   // method to compute first transport cross section per Volume (i.e. macroscropic first transport cross section; this 
159   // method is used only for testing and not d    160   // method is used only for testing and not during a normal simulation) 
160   G4double CrossSectionPerVolume(const G4Mater << 161   virtual G4double CrossSectionPerVolume(const G4Material*, const G4ParticleDefinition*, G4double kineticEnergy,
                                                   >> 162                                          G4double cutEnergy = 0.0, G4double maxEnergy = DBL_MAX);
161                                                   163 
162   void     StartTracking(G4Track*) override;   << 164   void     StartTracking(G4Track*);
163                                                   165 
164   void     SampleMSC();                           166   void     SampleMSC();
165                                                   167 
166   G4double GetTransportMeanFreePath(const G4Pa    168   G4double GetTransportMeanFreePath(const G4ParticleDefinition*, G4double);
167                                                   169 
168   void SetOptionPWACorrection(G4bool opt)    {    170   void SetOptionPWACorrection(G4bool opt)    { fIsUsePWACorrection = opt; }
169                                                   171 
170   G4bool GetOptionPWACorrection() const      {    172   G4bool GetOptionPWACorrection() const      { return fIsUsePWACorrection; }
171                                                   173 
172   void   SetOptionMottCorrection(G4bool opt) {    174   void   SetOptionMottCorrection(G4bool opt) { fIsUseMottCorrection = opt; }
173                                                   175 
174   G4bool GetOptionMottCorrection() const     {    176   G4bool GetOptionMottCorrection() const     { return fIsUseMottCorrection; }
175                                                   177 
176   G4GoudsmitSaundersonTable* GetGSTable()         178   G4GoudsmitSaundersonTable* GetGSTable()          { return fGSTable; }
177                                                   179 
178   G4GSPWACorrections*        GetPWACorrection(    180   G4GSPWACorrections*        GetPWACorrection()    { return fPWACorrection; }
179                                                   181 
180   //  hide assignment operator                 << 
181   G4GoudsmitSaundersonMscModel & operator=(con << 
182   G4GoudsmitSaundersonMscModel(const  G4Goudsm << 
183                                                << 
184 private:                                          182 private:
185   inline void     SetParticle(const G4Particle    183   inline void     SetParticle(const G4ParticleDefinition* p);
186                                                   184 
187   inline G4double GetLambda(G4double);            185   inline G4double GetLambda(G4double);
188                                                   186 
                                                   >> 187   //  hide assignment operator
                                                   >> 188   G4GoudsmitSaundersonMscModel & operator=(const  G4GoudsmitSaundersonMscModel &right);
                                                   >> 189   G4GoudsmitSaundersonMscModel(const  G4GoudsmitSaundersonMscModel&);
189   G4double GetTransportMeanFreePathOnly(const     190   G4double GetTransportMeanFreePathOnly(const G4ParticleDefinition*,G4double);
190                                                   191 
191   inline G4double Randomizetlimit();              192   inline G4double Randomizetlimit();
192                                                   193 
193 private:                                          194 private:
194   CLHEP::HepRandomEngine* rndmEngineMod;          195   CLHEP::HepRandomEngine* rndmEngineMod;
195   //                                              196   //
196   G4double currentKinEnergy;                      197   G4double currentKinEnergy;
197   G4double currentRange;                          198   G4double currentRange;
198   //                                              199   //
199   G4double fr;                                    200   G4double fr;
200   G4double rangeinit;                             201   G4double rangeinit;
201   G4double geombig;                               202   G4double geombig;
202   G4double geomlimit;                             203   G4double geomlimit;
203   G4double tlimit;                                204   G4double tlimit;
204   G4double tgeom;                                 205   G4double tgeom;
205   //                                              206   //
206   G4double par1;                                  207   G4double par1;
207   G4double par2;                                  208   G4double par2;
208   G4double par3;                                  209   G4double par3;
209   G4double tlimitminfix2;                         210   G4double tlimitminfix2;
210   G4double tausmall;                              211   G4double tausmall;
211   G4double mass;                                  212   G4double mass;
212   G4double taulim;                                213   G4double taulim;
213   //                                              214   //
214   //                                              215   //
215   G4double presafety;                             216   G4double presafety;
216   G4double fZeff;                                 217   G4double fZeff;
217   //                                              218   //
218   G4int    charge;                                219   G4int    charge;
219   G4int    currentMaterialIndex;                  220   G4int    currentMaterialIndex;
220   //                                              221   //
221   G4bool   firstStep;                             222   G4bool   firstStep;
222   //                                              223   //
223   G4LossTableManager*         theManager;         224   G4LossTableManager*         theManager;
224   const G4ParticleDefinition* particle;           225   const G4ParticleDefinition* particle;
225   G4ParticleChangeForMSC*     fParticleChange;    226   G4ParticleChangeForMSC*     fParticleChange;
226   const G4MaterialCutsCouple* currentCouple;      227   const G4MaterialCutsCouple* currentCouple;
227                                                   228 
228   G4GoudsmitSaundersonTable*  fGSTable;           229   G4GoudsmitSaundersonTable*  fGSTable;
229   G4GSPWACorrections*         fPWACorrection;     230   G4GSPWACorrections*         fPWACorrection;
230                                                   231 
231   G4bool   fIsUsePWACorrection;                   232   G4bool   fIsUsePWACorrection;
232   G4bool   fIsUseMottCorrection;                  233   G4bool   fIsUseMottCorrection;
233   //                                              234   //
234   G4double fLambda0; // elastic mean free path    235   G4double fLambda0; // elastic mean free path
235   G4double fLambda1; // first transport mean f    236   G4double fLambda1; // first transport mean free path
236   G4double fScrA;    // screening parameter       237   G4double fScrA;    // screening parameter
237   G4double fG1;      // first transport coef.     238   G4double fG1;      // first transport coef.
238   // in case of Mott-correction                   239   // in case of Mott-correction
239   G4double fMCtoScrA;                             240   G4double fMCtoScrA;
240   G4double fMCtoQ1;                               241   G4double fMCtoQ1;
241   G4double fMCtoG2PerG1;                          242   G4double fMCtoG2PerG1;
242   //                                              243   //
243   G4double fTheTrueStepLenght;                    244   G4double fTheTrueStepLenght;
244   G4double fTheTransportDistance;                 245   G4double fTheTransportDistance;
245   G4double fTheZPathLenght;                       246   G4double fTheZPathLenght;
246   //                                              247   //
247   G4ThreeVector fTheDisplacementVector;           248   G4ThreeVector fTheDisplacementVector;
248   G4ThreeVector fTheNewDirection;                 249   G4ThreeVector fTheNewDirection;
249   //                                              250   //
250   G4bool fIsEndedUpOnBoundary;  // step ended     251   G4bool fIsEndedUpOnBoundary;  // step ended up on boundary i.e. transportation is the winer
251   G4bool fIsMultipleSacettring;                   252   G4bool fIsMultipleSacettring;
252   G4bool fIsSingleScattering;                     253   G4bool fIsSingleScattering;
253   G4bool fIsEverythingWasDone;                    254   G4bool fIsEverythingWasDone;
254   G4bool fIsNoScatteringInMSC;                    255   G4bool fIsNoScatteringInMSC;
255   G4bool fIsNoDisplace;                           256   G4bool fIsNoDisplace;
256   G4bool fIsInsideSkin;                           257   G4bool fIsInsideSkin;
257   G4bool fIsWasOnBoundary;                        258   G4bool fIsWasOnBoundary;
258   G4bool fIsFirstRealStep;                        259   G4bool fIsFirstRealStep;
259   //                                              260   //
260   static G4bool gIsUseAccurate;                   261   static G4bool gIsUseAccurate;
261   static G4bool gIsOptimizationOn;                262   static G4bool gIsOptimizationOn;
262 };                                                263 };
263                                                   264 
264 //////////////////////////////////////////////    265 ////////////////////////////////////////////////////////////////////////////////
265 inline                                            266 inline
266 void G4GoudsmitSaundersonMscModel::SetParticle    267 void G4GoudsmitSaundersonMscModel::SetParticle(const G4ParticleDefinition* p)
267 {                                                 268 {
268   if (p != particle) {                            269   if (p != particle) {
269     particle = p;                                 270     particle = p;
270     charge = (G4int)(p->GetPDGCharge()/CLHEP::    271     charge = (G4int)(p->GetPDGCharge()/CLHEP::eplus);
271     mass = p->GetPDGMass();                       272     mass = p->GetPDGMass();
272   }                                               273   }
273 }                                                 274 }
274                                                   275 
275                                                   276 
276 //....oooOO0OOooo........oooOO0OOooo........oo    277 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
277 inline                                            278 inline
278 G4double G4GoudsmitSaundersonMscModel::Randomi    279 G4double G4GoudsmitSaundersonMscModel::Randomizetlimit()
279 {                                                 280 {
280   G4double temptlimit;                         << 281   G4double temptlimit = tlimit;
281     do {                                          282     do {
282          temptlimit = G4RandGauss::shoot(rndmE    283          temptlimit = G4RandGauss::shoot(rndmEngineMod,tlimit,0.1*tlimit);
283        } while ( (temptlimit<0.) || (temptlimi    284        } while ( (temptlimit<0.) || (temptlimit>2.*tlimit));
284                                                   285 
285   return temptlimit;                              286   return temptlimit;
286 }                                                 287 }
287                                                   288 
288                                                   289 
289                                                   290 
290 #endif                                            291 #endif
291                                                   292