Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/utils/src/G4VMultipleScattering.cc

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/utils/src/G4VMultipleScattering.cc (Version 11.3.0) and /processes/electromagnetic/utils/src/G4VMultipleScattering.cc (Version 6.0)


  1 //                                                  1 //
  2 // *******************************************      2 // ********************************************************************
  3 // * License and Disclaimer                    <<   3 // * DISCLAIMER                                                       *
  4 // *                                                4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of th <<   5 // * The following disclaimer summarizes all the specific disclaimers *
  6 // * the Geant4 Collaboration.  It is provided <<   6 // * of contributors to this software. The specific disclaimers,which *
  7 // * conditions of the Geant4 Software License <<   7 // * govern, are listed with their locations in:                      *
  8 // * LICENSE and available at  http://cern.ch/ <<   8 // *   http://cern.ch/geant4/license                                  *
  9 // * include a list of copyright holders.      << 
 10 // *                                                9 // *                                                                  *
 11 // * Neither the authors of this software syst     10 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing fin     11 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warran     12 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assum     13 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  <<  14 // * use.                                                             *
 16 // * for the full disclaimer and the limitatio << 
 17 // *                                               15 // *                                                                  *
 18 // * This  code  implementation is the result  <<  16 // * This  code  implementation is the  intellectual property  of the *
 19 // * technical work of the GEANT4 collaboratio <<  17 // * GEANT4 collaboration.                                            *
 20 // * By using,  copying,  modifying or  distri <<  18 // * By copying,  distributing  or modifying the Program (or any work *
 21 // * any work based  on the software)  you  ag <<  19 // * based  on  the Program)  you indicate  your  acceptance of  this *
 22 // * use  in  resulting  scientific  publicati <<  20 // * statement, and all its terms.                                    *
 23 // * acceptance of all terms of the Geant4 Sof << 
 24 // *******************************************     21 // ********************************************************************
 25 //                                                 22 //
                                                   >>  23 // $Id: G4VMultipleScattering.cc,v 1.21 2003/11/11 14:05:10 vnivanch Exp $
                                                   >>  24 // GEANT4 tag $Name: geant4-06-00 $
 26 //                                                 25 //
 27 // -------------------------------------------     26 // -------------------------------------------------------------------
 28 //                                                 27 //
 29 // GEANT4 Class file                               28 // GEANT4 Class file
 30 //                                                 29 //
 31 //                                                 30 //
 32 // File name:     G4VMultipleScattering            31 // File name:     G4VMultipleScattering
 33 //                                                 32 //
 34 // Author:        Vladimir Ivanchenko on base      33 // Author:        Vladimir Ivanchenko on base of Laszlo Urban code
 35 //                                                 34 //
 36 // Creation date: 25.03.2003                       35 // Creation date: 25.03.2003
 37 //                                                 36 //
 38 // Modifications:                                  37 // Modifications:
 39 //                                                 38 //
                                                   >>  39 // 13.04.03 Change printout (V.Ivanchenko)
                                                   >>  40 // 04-06-03 Fix compilation warnings (V.Ivanchenko)
 40 // 16-07-03 Use G4VMscModel interface (V.Ivanc     41 // 16-07-03 Use G4VMscModel interface (V.Ivanchenko)
 41 // 03-11-03 Fix initialisation problem in Retr     42 // 03-11-03 Fix initialisation problem in RetrievePhysicsTable (V.Ivanchenko)
 42 // 04-11-03 Update PrintInfoDefinition (V.Ivan     43 // 04-11-03 Update PrintInfoDefinition (V.Ivanchenko)
 43 // 01-03-04 SampleCosineTheta signature change << 
 44 // 22-04-04 SampleCosineTheta signature change << 
 45 // 27-08-04 Add InitialiseForRun method (V.Iva << 
 46 // 08-11-04 Migration to new interface of Stor << 
 47 // 11-03-05 Shift verbose level by 1 (V.Ivantc << 
 48 // 15-04-05 optimize internal interface (V.Iva << 
 49 // 15-04-05 remove boundary flag (V.Ivanchenko << 
 50 // 27-10-05 introduce virtual function MscStep << 
 51 // 12-04-07 Add verbosity at destruction (V.Iv << 
 52 // 27-10-07 Virtual functions moved to source  << 
 53 // 11-03-08 Set skin value does not effect ste << 
 54 // 24-06-09 Removed hidden bin in G4PhysicsVec << 
 55 // 04-06-13 Adoptation to MT mode (V.Ivanchenk << 
 56 //                                                 44 //
                                                   >>  45 //
                                                   >>  46 // Class Description:
                                                   >>  47 //
                                                   >>  48 // It is the generic process of multiple scattering it includes common
                                                   >>  49 // part of calculations for all charged particles
 57                                                    50 
 58 // -------------------------------------------     51 // -------------------------------------------------------------------
 59 //                                                 52 //
 60 //....oooOO0OOooo........oooOO0OOooo........oo     53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 61 //....oooOO0OOooo........oooOO0OOooo........oo     54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 62                                                    55 
 63 #include "G4VMultipleScattering.hh"                56 #include "G4VMultipleScattering.hh"
 64 #include "G4PhysicalConstants.hh"              << 
 65 #include "G4SystemOfUnits.hh"                  << 
 66 #include "G4LossTableManager.hh"                   57 #include "G4LossTableManager.hh"
 67 #include "G4MaterialCutsCouple.hh"             << 
 68 #include "G4Step.hh"                               58 #include "G4Step.hh"
 69 #include "G4ParticleDefinition.hh"                 59 #include "G4ParticleDefinition.hh"
 70 #include "G4VEmFluctuationModel.hh"                60 #include "G4VEmFluctuationModel.hh"
                                                   >>  61 #include "G4DataVector.hh"
                                                   >>  62 #include "G4PhysicsTable.hh"
                                                   >>  63 #include "G4PhysicsVector.hh"
                                                   >>  64 #include "G4PhysicsLogVector.hh"
 71 #include "G4UnitsTable.hh"                         65 #include "G4UnitsTable.hh"
 72 #include "G4ProductionCutsTable.hh"                66 #include "G4ProductionCutsTable.hh"
 73 #include "G4Electron.hh"                       <<  67 #include "G4Region.hh"
 74 #include "G4GenericIon.hh"                     <<  68 #include "G4RegionStore.hh"
                                                   >>  69 #include "G4Navigator.hh"
 75 #include "G4TransportationManager.hh"              70 #include "G4TransportationManager.hh"
 76 #include "G4SafetyHelper.hh"                   << 
 77 #include "G4ParticleTable.hh"                  << 
 78 #include "G4ProcessVector.hh"                  << 
 79 #include "G4ProcessManager.hh"                 << 
 80 #include "G4LossTableBuilder.hh"               << 
 81 #include "G4EmTableUtil.hh"                    << 
 82 #include <iostream>                            << 
 83                                                    71 
 84 //....oooOO0OOooo........oooOO0OOooo........oo     72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 85                                                    73 
 86 G4VMultipleScattering::G4VMultipleScattering(c <<  74 G4VMultipleScattering::G4VMultipleScattering(const G4String& name, G4ProcessType type):
 87   : G4VContinuousDiscreteProcess("msc", fElect <<  75                  G4VContinuousDiscreteProcess(name, type),
 88   fNewPosition(0.,0.,0.),                      <<  76   navigator(0),
 89   fNewDirection(0.,0.,1.)                      <<  77   theLambdaTable(0),
 90 {                                              <<  78   currentCouple(0),
 91   theParameters = G4EmParameters::Instance();  <<  79   nBins(120),
 92   SetVerboseLevel(1);                          <<  80   boundary(false),
 93   SetProcessSubType(fMultipleScattering);      <<  81   latDisplasment(true),
 94                                                <<  82   buildLambdaTable(true)
 95   lowestKinEnergy = 10*CLHEP::eV;              <<  83 {
 96                                                <<  84   minKinEnergy = 100.0*eV;
 97   geomMin = 0.05*CLHEP::nm;                    <<  85   maxKinEnergy = 100.0*TeV;
 98   minDisplacement2 = geomMin*geomMin;          <<  86   SetVerboseLevel(0);
 99                                                << 
100   pParticleChange = &fParticleChange;          << 
101                                                << 
102   modelManager = new G4EmModelManager();           87   modelManager = new G4EmModelManager();
103   emManager = G4LossTableManager::Instance();  <<  88   (G4LossTableManager::Instance())->Register(this);
104   mscModels.reserve(2);                        <<  89 
105   emManager->Register(this);                   << 
106 }                                                  90 }
107                                                    91 
108 //....oooOO0OOooo........oooOO0OOooo........oo     92 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
109                                                    93 
110 G4VMultipleScattering::~G4VMultipleScattering(     94 G4VMultipleScattering::~G4VMultipleScattering()
111 {                                                  95 {
                                                   >>  96   (G4LossTableManager::Instance())->DeRegister(this);
112   delete modelManager;                             97   delete modelManager;
113   emManager->DeRegister(this);                 <<  98   if (theLambdaTable) {
114 }                                              <<  99     theLambdaTable->clearAndDestroy();
115                                                << 100     delete theLambdaTable;
116 //....oooOO0OOooo........oooOO0OOooo........oo << 
117                                                << 
118 void G4VMultipleScattering::AddEmModel(G4int o << 
119                                        const G << 
120 {                                              << 
121   if(nullptr == ptr) { return; }               << 
122   G4VEmFluctuationModel* fm = nullptr;         << 
123   modelManager->AddEmModel(order, ptr, fm, reg << 
124   ptr->SetParticleChange(pParticleChange);     << 
125 }                                              << 
126                                                << 
127 //....oooOO0OOooo........oooOO0OOooo........oo << 
128                                                << 
129 void G4VMultipleScattering::SetEmModel(G4VMscM << 
130 {                                              << 
131   if(nullptr == ptr) { return; }               << 
132   if(!mscModels.empty()) {                     << 
133     for(auto & msc : mscModels) { if(msc == pt << 
134   }                                               101   }
135   mscModels.push_back(ptr);                    << 102   (G4LossTableManager::Instance())->DeRegister(this);
136 }                                                 103 }
137                                                   104 
138 //....oooOO0OOooo........oooOO0OOooo........oo    105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
139                                                   106 
140 void                                           << 107 void G4VMultipleScattering::BuildPhysicsTable(const G4ParticleDefinition& part)
141 G4VMultipleScattering::PreparePhysicsTable(con << 
142 {                                                 108 {
143   G4bool master = emManager->IsMaster();       << 109   currentParticle = &part;
144   if (nullptr == firstParticle) { firstParticl << 110   currentCouple = 0;
                                                   >> 111   if(0 < verboseLevel) {
                                                   >> 112     G4cout << "========================================================" << G4endl;
                                                   >> 113     G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
                                                   >> 114            << GetProcessName()
                                                   >> 115            << " and particle " << part.GetParticleName()
                                                   >> 116            << G4endl;
                                                   >> 117   }
                                                   >> 118 
                                                   >> 119   InitialiseProcess(part);
                                                   >> 120 
                                                   >> 121   if(latDisplasment) navigator = G4TransportationManager::GetTransportationManager()
                                                   >> 122              ->GetNavigatorForTracking();
                                                   >> 123   const G4DataVector* theCuts = modelManager->Initialise(&part, 0, 10.0, verboseLevel);
                                                   >> 124 
                                                   >> 125   if (buildLambdaTable) {
                                                   >> 126 
                                                   >> 127     const G4ProductionCutsTable* theCoupleTable=
                                                   >> 128           G4ProductionCutsTable::GetProductionCutsTable();
                                                   >> 129     size_t numOfCouples = theCoupleTable->GetTableSize();
                                                   >> 130     theLambdaTable = new G4PhysicsTable(numOfCouples);
                                                   >> 131 
                                                   >> 132     for (size_t i=0; i<numOfCouples; i++) {
                                                   >> 133 
                                                   >> 134       // create physics vector and fill it
                                                   >> 135       const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i);
                                                   >> 136       G4PhysicsVector* aVector = PhysicsVector(couple);
                                                   >> 137       modelManager->FillLambdaVector(aVector, couple, false);
145                                                   138 
146   emManager->PreparePhysicsTable(&part, this); << 139       // Insert vector for this material into the table
147   currParticle = nullptr;                      << 140       theLambdaTable->insert(aVector) ;
                                                   >> 141     }
148                                                   142 
149   if(firstParticle == &part) {                 << 143     if(0 < verboseLevel) {
150     baseMat = emManager->GetTableBuilder()->Ge << 144       G4cout << "Lambda table is built for "
151     G4EmTableUtil::PrepareMscProcess(this, par << 145              << part.GetParticleName()
152              stepLimit, facrange,              << 146              << G4endl;
153              latDisplacement, master,          << 
154              isIon, baseMat);                  << 
155                                                << 
156     numberOfModels = modelManager->NumberOfMod << 
157     currentModel = GetModelByIndex(0);         << 
158                                                << 
159     if (nullptr == safetyHelper) {             << 
160       safetyHelper = G4TransportationManager:: << 
161   ->GetSafetyHelper();                         << 
162       safetyHelper->InitialiseHelper();        << 
163     }                                             147     }
164   }                                            << 148     if(2 < verboseLevel) G4cout << *theLambdaTable << G4endl;
165 }                                              << 149     if(5 < verboseLevel) G4cout << theCuts << G4endl;
166                                                   150 
167 //....oooOO0OOooo........oooOO0OOooo........oo << 151   }
168                                                   152 
169 void G4VMultipleScattering::BuildPhysicsTable( << 153   G4String num = part.GetParticleName();
170 {                                              << 154   if (verboseLevel>0 || num == "e-" || num == "mu+" || num == "proton"
171   G4bool master = emManager->IsMaster();       << 155                      || num == "pi-") PrintInfoDefinition();
172                                                   156 
173   if(firstParticle == &part) {                 << 157   if(0 < verboseLevel) {
174     emManager->BuildPhysicsTable(&part);       << 158     G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
175   }                                            << 159            << GetProcessName()
176   const G4VMultipleScattering* ptr = this;     << 160            << " and particle " << part.GetParticleName()
177   if(!master) {                                << 161            << G4endl;
178     ptr = static_cast<const G4VMultipleScatter << 
179   }                                               162   }
180                                                << 
181   G4EmTableUtil::BuildMscProcess(this, ptr, pa << 
182          numberOfModels, master);              << 
183 }                                                 163 }
184                                                   164 
185 //....oooOO0OOooo........oooOO0OOooo........oo    165 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
186                                                   166 
187 void G4VMultipleScattering::StreamInfo(std::os << 167 void G4VMultipleScattering::AddEmModel(G4int order, G4VEmModel* p,
188                   const G4ParticleDefinition&  << 168                                  const G4Region* region)
189 {                                                 169 {
190   G4String indent = (rst ? "  " : "");         << 170   G4VEmFluctuationModel* fm = 0;
191   outFile << G4endl << indent << GetProcessNam << 171   modelManager->AddEmModel(order, p, fm, region);
192   if (!rst) outFile << " for " << part.GetPart << 
193   outFile  << "  SubType= " << GetProcessSubTy << 
194   modelManager->DumpModelList(outFile, verbose << 
195 }                                                 172 }
196                                                   173 
197 //....oooOO0OOooo........oooOO0OOooo........oo    174 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
198                                                   175 
199 void G4VMultipleScattering::StartTracking(G4Tr << 176 G4VParticleChange* G4VMultipleScattering::PostStepDoIt(const G4Track& track,
                                                   >> 177                                                        const G4Step& step)
200 {                                                 178 {
201   G4VEnergyLossProcess* eloss = nullptr;       << 179   fParticleChange.Initialize(track);
202   if(track->GetParticleDefinition() != currPar << 180   G4double kineticEnergy = track.GetKineticEnergy();
203     currParticle = track->GetParticleDefinitio << 181   G4double truestep = step.GetStepLength();
204     fIonisation = emManager->GetEnergyLossProc << 
205     eloss = fIonisation;                       << 
206   }                                            << 
207   for(G4int i=0; i<numberOfModels; ++i) {      << 
208     G4VMscModel* msc = GetModelByIndex(i);     << 
209     msc->StartTracking(track);                 << 
210     if(nullptr != eloss) {                     << 
211       msc->SetIonisation(eloss, currParticle); << 
212     }                                          << 
213   }                                            << 
214 }                                              << 
215                                                   182 
216 //....oooOO0OOooo........oooOO0OOooo........oo << 183   if (kineticEnergy > 0.0) {
                                                   >> 184     G4double cth  = currentModel->SampleCosineTheta(truestep);
                                                   >> 185     G4double sth  = sqrt(1.-cth*cth);
                                                   >> 186     G4double phi  = twopi*G4UniformRand();
                                                   >> 187     G4double dirx = sth*cos(phi);
                                                   >> 188     G4double diry = sth*sin(phi);
                                                   >> 189 
                                                   >> 190     G4ThreeVector oldDirection = track.GetMomentumDirection();
                                                   >> 191     G4ThreeVector newDirection(dirx,diry,cth);
                                                   >> 192     newDirection.rotateUz(oldDirection);
                                                   >> 193     fParticleChange.SetMomentumChange(newDirection);
217                                                   194 
218 G4double G4VMultipleScattering::AlongStepGetPh << 
219                              const G4Track& tr << 
220                              G4double,         << 
221                              G4double currentM << 
222                              G4double&,        << 
223                              G4GPILSelection*  << 
224 {                                              << 
225   // get Step limit proposed by the process    << 
226   *selection = NotCandidateForSelection;       << 
227   physStepLimit = gPathLength = tPathLength =  << 
228                                                   195 
229   G4double ekin = track.GetKineticEnergy();    << 
230   /*                                              196   /*
231   G4cout << "MSC::AlongStepGPIL: Ekin= " << ek << 197   if(0 < verboseLevel) {
232          << "  " << currParticle->GetParticleN << 198     const G4ParticleDefinition* pd = dynParticle->GetDefinition();
233          << " currMod " << currentModel        << 199     G4cout << "G4VMultipleScattering::PostStepDoIt: Sample secondary; E= " << finalT/MeV
234          << G4endl;                            << 200            << " MeV; model= (" << currentModel->LowEnergyLimit(pd)
235   */                                           << 201            << ", " <<  currentModel->HighEnergyLimit(pd) << ")"
236   // isIon flag is used only to select a model << 202            << G4endl;
237   if(isIon) {                                  << 
238     ekin *= proton_mass_c2/track.GetParticleDe << 
239   }                                            << 
240   const G4MaterialCutsCouple* couple = track.G << 
241                                                << 
242   // select new model, static cast is possible << 
243   if(1 < numberOfModels) {                     << 
244     currentModel =                             << 
245       static_cast<G4VMscModel*>(SelectModel(ek << 
246   }                                            << 
247   currentModel->SetCurrentCouple(couple);      << 
248   // msc is active is model is active, energy  << 
249   // and step size is above the limit;         << 
250   // if it is active msc may limit the step    << 
251   if(currentModel->IsActive(ekin) && tPathLeng << 
252      && ekin >= lowestKinEnergy) {             << 
253     isActive = true;                           << 
254     tPathLength =                              << 
255       currentModel->ComputeTruePathLengthLimit << 
256     if (tPathLength < physStepLimit) {         << 
257       *selection = CandidateForSelection;      << 
258     }                                          << 
259   } else {                                     << 
260     isActive = false;                          << 
261     gPathLength = DBL_MAX;                     << 
262   }                                               203   }
263                                                << 
264   //if(currParticle->GetPDGMass() > GeV)       << 
265   /*                                           << 
266   G4cout << "MSC::AlongStepGPIL: Ekin= " << ek << 
267          << "  " << currParticle->GetParticleN << 
268          << " gPathLength= " << gPathLength    << 
269          << " tPathLength= " << tPathLength    << 
270          << " currentMinimalStep= " << current << 
271          << " isActive " << isActive << G4endl << 
272   */                                              204   */
273   return gPathLength;                          << 
274 }                                              << 
275                                                   205 
276 //....oooOO0OOooo........oooOO0OOooo........oo << 206     //    G4cout << "PostStep: sth= " << sth << " trueLength= " << truestep << " tLast= " << truePathLength << G4endl;
277                                                   207 
278 G4double                                       << 208     if (latDisplasment) {
279 G4VMultipleScattering::PostStepGetPhysicalInte << 
280               const G4Track&, G4double, G4Forc << 
281 {                                              << 
282   *condition = NotForced;                      << 
283   return DBL_MAX;                              << 
284 }                                              << 
285                                                   209 
286 //....oooOO0OOooo........oooOO0OOooo........oo << 210       G4double safety = step.GetPostStepPoint()->GetSafety();
                                                   >> 211       if ( safety > 0.0) {
                                                   >> 212         G4double r = currentModel->SampleDisplacement();
                                                   >> 213         if (r > safety) r = safety;
287                                                   214 
288 G4VParticleChange*                             << 215   //    G4cout << "r= " << r << " safety= " << safety << G4endl;
289 G4VMultipleScattering::AlongStepDoIt(const G4T << 
290 {                                              << 
291   fParticleChange.InitialiseMSC(track, step);  << 
292   fNewPosition = fParticleChange.GetProposedPo << 
293   fPositionChanged = false;                    << 
294                                                   216 
295   G4double geomLength = step.GetStepLength();  << 217         // sample direction of lateral displacement
                                                   >> 218         G4double phi  = twopi*G4UniformRand();
                                                   >> 219         G4double dirx = cos(phi);
                                                   >> 220         G4double diry = sin(phi);
296                                                   221 
297   // very small step - no msc                  << 222         G4ThreeVector latDirection(dirx,diry,0.0);
298   if(!isActive) {                              << 223         latDirection.rotateUz(oldDirection);
299     tPathLength = geomLength;                  << 
300                                                   224 
301     // sample msc                              << 225         // compute new endpoint of the Step
302   } else {                                     << 226         G4ThreeVector newPosition = (step.GetPostStepPoint())->GetPosition()
303     G4double range =                           << 227               + r*latDirection;
304       currentModel->GetRange(currParticle,trac << 228 
305                              track.GetMaterial << 229         navigator->LocateGlobalPointWithinVolume(newPosition);
306                                                << 230 
307     tPathLength = currentModel->ComputeTrueSte << 231         fParticleChange.SetPositionChange(newPosition);
308                                                << 
309     /*                                         << 
310     if(currParticle->GetPDGMass() > 0.9*GeV)   << 
311     G4cout << "G4VMsc::AlongStepDoIt: GeomLeng << 
312            << geomLength                       << 
313            << " tPathLength= " << tPathLength  << 
314            << " physStepLimit= " << physStepLi << 
315            << " dr= " << range - tPathLength   << 
316            << " ekin= " << track.GetKineticEne << 
317     */                                         << 
318     // protection against wrong t->g->t conver << 
319     tPathLength = std::min(tPathLength, physSt << 
320                                                << 
321     // do not sample scattering at the last or << 
322     if(tPathLength < range && tPathLength > ge << 
323                                                << 
324       static const G4double minSafety = 1.20*C << 
325       static const G4double sFact = 0.99;      << 
326                                                << 
327       G4ThreeVector displacement = currentMode << 
328         step.GetPostStepPoint()->GetMomentumDi << 
329                                                << 
330       G4double r2 = displacement.mag2();       << 
331       //G4cout << "    R= " << sqrt(r2) << " R << 
332       //     << " flag= " << fDispBeyondSafety << 
333       if(r2 > minDisplacement2) {              << 
334                                                << 
335         fPositionChanged = true;               << 
336         G4double dispR = std::sqrt(r2);        << 
337         G4double postSafety =                  << 
338           sFact*safetyHelper->ComputeSafety(fN << 
339         //G4cout<<"    R= "<< dispR<<" postSaf << 
340                                                << 
341         // far away from geometry boundary     << 
342         if(postSafety > 0.0 && dispR <= postSa << 
343           fNewPosition += displacement;        << 
344                                                << 
345           //near the boundary                  << 
346         } else {                               << 
347           // displaced point is definitely wit << 
348           //G4cout<<"    R= "<<dispR<<" postSa << 
349           if(dispR < postSafety) {             << 
350             fNewPosition += displacement;      << 
351                                                << 
352             // reduced displacement            << 
353           } else if(postSafety > geomMin) {    << 
354             fNewPosition += displacement*(post << 
355                                                << 
356             // very small postSafety           << 
357           } else {                             << 
358             fPositionChanged = false;          << 
359           }                                    << 
360         }                                      << 
361         if(fPositionChanged) {                 << 
362           safetyHelper->ReLocateWithinVolume(f << 
363           fParticleChange.ProposePosition(fNew << 
364         }                                      << 
365       }                                           232       }
366     }                                             233     }
367   }                                               234   }
368   fParticleChange.ProposeTrueStepLength(tPathL << 
369   return &fParticleChange;                     << 
370 }                                              << 
371                                                   235 
372 //....oooOO0OOooo........oooOO0OOooo........oo << 236   return &fParticleChange;
373                                                << 
374 G4double G4VMultipleScattering::GetContinuousS << 
375                                        const G << 
376                                        G4doubl << 
377                                        G4doubl << 
378                                        G4doubl << 
379 {                                              << 
380   G4GPILSelection selection = NotCandidateForS << 
381   G4double x = AlongStepGetPhysicalInteraction << 
382                                                << 
383                                                << 
384                                                << 
385   return x;                                    << 
386 }                                                 237 }
387                                                   238 
388 //....oooOO0OOooo........oooOO0OOooo........oo    239 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
389                                                   240 
390 G4double G4VMultipleScattering::ContinuousStep << 241 void G4VMultipleScattering::PrintInfoDefinition()
391                                        const G << 
392                                        G4doubl << 
393                                        G4doubl << 
394                                        G4doubl << 
395 {                                                 242 {
396   return GetContinuousStepLimit(track,previous << 243   G4cout << G4endl << GetProcessName() << ":  Model variant of multiple scattering " 
397                                 currentSafety) << 244          << "for " << currentParticle->GetParticleName()
                                                   >> 245          << G4endl;
                                                   >> 246   if (theLambdaTable) {
                                                   >> 247     G4cout << "      Lambda tables from "
                                                   >> 248            << G4BestUnit(MinKinEnergy(),"Energy")
                                                   >> 249      << " to "
                                                   >> 250            << G4BestUnit(MaxKinEnergy(),"Energy")
                                                   >> 251            << " in " << nBins << " bins."
                                                   >> 252            << G4endl;
                                                   >> 253   }
                                                   >> 254   if (1 < verboseLevel) {
                                                   >> 255       G4cout << "LambdaTable address= " << theLambdaTable << G4endl;
                                                   >> 256       if(theLambdaTable) G4cout << (*theLambdaTable) << G4endl;
                                                   >> 257   }
398 }                                                 258 }
399                                                   259 
400 //....oooOO0OOooo........oooOO0OOooo........oo    260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
401                                                   261 
402 G4double G4VMultipleScattering::GetMeanFreePat << 262 G4PhysicsVector* G4VMultipleScattering::PhysicsVector(const G4MaterialCutsCouple* couple)
403               const G4Track&, G4double, G4Forc << 
404 {                                                 263 {
405   *condition = Forced;                         << 264   G4int nbins = 3;
406   return DBL_MAX;                              << 265   //G4int nbins =  nDEDXBins;
                                                   >> 266   if( couple->IsUsed() ) nbins = nBins;
                                                   >> 267   //  G4double xmax = maxKinEnergy*exp( log(maxKinEnergy/minKinEnergy) / ((G4double)(nbins-1)) );
                                                   >> 268   G4PhysicsVector* v = new G4PhysicsLogVector(minKinEnergy, maxKinEnergy, nbins);
                                                   >> 269   return v;
407 }                                                 270 }
408                                                   271 
409 //....oooOO0OOooo........oooOO0OOooo........oo << 272 G4bool G4VMultipleScattering::StorePhysicsTable(G4ParticleDefinition* part,
410                                                << 273              const G4String& directory,
411 G4bool                                         << 274                    G4bool ascii)
412 G4VMultipleScattering::StorePhysicsTable(const << 
413                                          const << 
414                                          G4boo << 
415 {                                                 275 {
416   G4bool yes = true;                           << 276   G4bool res = true;
417   if(part != firstParticle || !emManager->IsMa << 277   if ( theLambdaTable ) {
                                                   >> 278     const G4String name = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
                                                   >> 279     G4bool yes = theLambdaTable->StorePhysicsTable(name,ascii);
                                                   >> 280     if( !yes ) res = false;
                                                   >> 281   }
418                                                   282 
419   return G4EmTableUtil::StoreMscTable(this, pa << 283   if ( res ) {
420               numberOfModels, verboseLevel,    << 284     G4cout << "Physics table are stored for " << part->GetParticleName()
421                                       ascii);  << 285            << " and process " << GetProcessName()
                                                   >> 286      << " in the directory <" << directory
                                                   >> 287      << "> " << G4endl;
                                                   >> 288   } else {
                                                   >> 289     G4cout << "Fail to store Physics Table for " << part->GetParticleName()
                                                   >> 290            << " and process " << GetProcessName()
                                                   >> 291      << " in the directory <" << directory
                                                   >> 292      << "> " << G4endl;
                                                   >> 293   }
                                                   >> 294   return res;
422 }                                                 295 }
423                                                   296 
424 //....oooOO0OOooo........oooOO0OOooo........oo    297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
425                                                   298 
426 G4bool                                         << 299 G4bool G4VMultipleScattering::RetrievePhysicsTable(G4ParticleDefinition* part,
427 G4VMultipleScattering::RetrievePhysicsTable(co << 300                const G4String& directory,
428                                             co << 301                              G4bool ascii)
429                                             G4 << 302 {
430 {                                              << 303   currentParticle = part;
431   return true;                                 << 304   if(0 < verboseLevel) {
432 }                                              << 305     G4cout << "========================================================" << G4endl;
433                                                << 306     G4cout << "G4VMultipleScattering::RetrievePhysicsTable() for "
434 //....oooOO0OOooo........oooOO0OOooo........oo << 307            << part->GetParticleName() << " and process "
435                                                << 308      << GetProcessName() << G4endl;
436 void G4VMultipleScattering::ProcessDescription << 309   }
437 {                                              << 310 
438   if(nullptr != firstParticle) {               << 311   InitialiseProcess(*part);
439     StreamInfo(outFile, *firstParticle, true); << 312   if(latDisplasment) navigator = G4TransportationManager::GetTransportationManager()
                                                   >> 313                                ->GetNavigatorForTracking();
                                                   >> 314 
                                                   >> 315   const G4DataVector* theCuts = modelManager->Initialise(part, 0, 10.0, verboseLevel);
                                                   >> 316   if(5 < verboseLevel) G4cout << theCuts << G4endl;
                                                   >> 317   if(!buildLambdaTable) return true;
                                                   >> 318 
                                                   >> 319   G4String num = part->GetParticleName();
                                                   >> 320   const G4ProductionCutsTable* theCoupleTable=
                                                   >> 321           G4ProductionCutsTable::GetProductionCutsTable();
                                                   >> 322   size_t numOfCouples = theCoupleTable->GetTableSize();
                                                   >> 323 
                                                   >> 324   G4String filename = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
                                                   >> 325   theLambdaTable = new G4PhysicsTable(numOfCouples);
                                                   >> 326   G4bool res = theLambdaTable->RetrievePhysicsTable(filename,ascii);
                                                   >> 327   if ( res ) {
                                                   >> 328     if (-1 < verboseLevel) {
                                                   >> 329         G4cout << "Lambda table for " << num << " is retrieved from <"
                                                   >> 330                << filename << ">"
                                                   >> 331                << G4endl;
                                                   >> 332     }
                                                   >> 333   } else {
                                                   >> 334     theLambdaTable->clearAndDestroy();
                                                   >> 335     theLambdaTable = 0;
                                                   >> 336     if (-1 < verboseLevel) {
                                                   >> 337         G4cout << "Lambda table for " << num << " in file <"
                                                   >> 338                << filename << "> is not exist"
                                                   >> 339                << G4endl;
                                                   >> 340     }
440   }                                               341   }
                                                   >> 342 
                                                   >> 343   if (verboseLevel>0 || num == "e-" || num == "mu+" || num == "proton")
                                                   >> 344            PrintInfoDefinition();
                                                   >> 345   return res;
441 }                                                 346 }
442                                                   347 
443 //....oooOO0OOooo........oooOO0OOooo........oo    348 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
444                                                   349 
445                                                   350