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 5.1.p1)


  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 //
 26 //                                                 23 //
 27 // -------------------------------------------     24 // -------------------------------------------------------------------
 28 //                                                 25 //
 29 // GEANT4 Class file                               26 // GEANT4 Class file
 30 //                                                 27 //
 31 //                                                 28 //
 32 // File name:     G4VMultipleScattering            29 // File name:     G4VMultipleScattering
 33 //                                                 30 //
 34 // Author:        Vladimir Ivanchenko on base      31 // Author:        Vladimir Ivanchenko on base of Laszlo Urban code
 35 //                                                 32 //
 36 // Creation date: 25.03.2003                       33 // Creation date: 25.03.2003
 37 //                                                 34 //
 38 // Modifications:                              <<  35 // Modifications: 13.04.03 Change printout (V.Ivanchenko)
 39 //                                                 36 //
 40 // 16-07-03 Use G4VMscModel interface (V.Ivanc << 
 41 // 03-11-03 Fix initialisation problem in Retr << 
 42 // 04-11-03 Update PrintInfoDefinition (V.Ivan << 
 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 //                                                 37 //
                                                   >>  38 // Class Description:
                                                   >>  39 //
                                                   >>  40 // It is the generic process of multiple scattering it includes common
                                                   >>  41 // part of calculations for all charged particles
 57                                                    42 
 58 // -------------------------------------------     43 // -------------------------------------------------------------------
 59 //                                                 44 //
 60 //....oooOO0OOooo........oooOO0OOooo........oo     45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 61 //....oooOO0OOooo........oooOO0OOooo........oo     46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 62                                                    47 
 63 #include "G4VMultipleScattering.hh"                48 #include "G4VMultipleScattering.hh"
 64 #include "G4PhysicalConstants.hh"              << 
 65 #include "G4SystemOfUnits.hh"                  << 
 66 #include "G4LossTableManager.hh"                   49 #include "G4LossTableManager.hh"
 67 #include "G4MaterialCutsCouple.hh"             << 
 68 #include "G4Step.hh"                               50 #include "G4Step.hh"
 69 #include "G4ParticleDefinition.hh"                 51 #include "G4ParticleDefinition.hh"
                                                   >>  52 #include "G4VEmModel.hh"
 70 #include "G4VEmFluctuationModel.hh"                53 #include "G4VEmFluctuationModel.hh"
                                                   >>  54 #include "G4DataVector.hh"
                                                   >>  55 #include "G4PhysicsTable.hh"
                                                   >>  56 #include "G4PhysicsVector.hh"
                                                   >>  57 #include "G4PhysicsLogVector.hh"
 71 #include "G4UnitsTable.hh"                         58 #include "G4UnitsTable.hh"
 72 #include "G4ProductionCutsTable.hh"                59 #include "G4ProductionCutsTable.hh"
 73 #include "G4Electron.hh"                       <<  60 #include "G4Region.hh"
 74 #include "G4GenericIon.hh"                     <<  61 #include "G4RegionStore.hh"
                                                   >>  62 #include "G4Navigator.hh"
 75 #include "G4TransportationManager.hh"              63 #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                                                    64 
 84 //....oooOO0OOooo........oooOO0OOooo........oo     65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 85                                                    66 
 86 G4VMultipleScattering::G4VMultipleScattering(c <<  67 G4VMultipleScattering::G4VMultipleScattering(const G4String& name, G4ProcessType type):
 87   : G4VContinuousDiscreteProcess("msc", fElect <<  68                  G4VContinuousDiscreteProcess(name, type),
 88   fNewPosition(0.,0.,0.),                      <<  69   navigator(0),
 89   fNewDirection(0.,0.,1.)                      <<  70   theLambdaTable(0),
                                                   >>  71   currentCouple(0),
                                                   >>  72   nBins(110),
                                                   >>  73   minKinEnergy(1.0*eV),
                                                   >>  74   maxKinEnergy(100.0*GeV),
                                                   >>  75   boundary(false),
                                                   >>  76   latDisplasment(true),
                                                   >>  77   buildLambdaTable(true)
 90 {                                                  78 {
 91   theParameters = G4EmParameters::Instance();  << 
 92   SetVerboseLevel(1);                          << 
 93   SetProcessSubType(fMultipleScattering);      << 
 94                                                << 
 95   lowestKinEnergy = 10*CLHEP::eV;              << 
 96                                                << 
 97   geomMin = 0.05*CLHEP::nm;                    << 
 98   minDisplacement2 = geomMin*geomMin;          << 
 99                                                << 
100   pParticleChange = &fParticleChange;          << 
101                                                << 
102   modelManager = new G4EmModelManager();           79   modelManager = new G4EmModelManager();
103   emManager = G4LossTableManager::Instance();  <<  80   (G4LossTableManager::Instance())->Register(this);
104   mscModels.reserve(2);                        <<  81 
105   emManager->Register(this);                   << 
106 }                                                  82 }
107                                                    83 
108 //....oooOO0OOooo........oooOO0OOooo........oo     84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
109                                                    85 
110 G4VMultipleScattering::~G4VMultipleScattering(     86 G4VMultipleScattering::~G4VMultipleScattering()
111 {                                                  87 {
                                                   >>  88   (G4LossTableManager::Instance())->DeRegister(this);
112   delete modelManager;                             89   delete modelManager;
113   emManager->DeRegister(this);                 <<  90   if (theLambdaTable) {
                                                   >>  91     theLambdaTable->clearAndDestroy();
                                                   >>  92     delete theLambdaTable;
                                                   >>  93   }
                                                   >>  94   (G4LossTableManager::Instance())->DeRegister(this);
114 }                                                  95 }
115                                                    96 
116 //....oooOO0OOooo........oooOO0OOooo........oo     97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
117                                                    98 
118 void G4VMultipleScattering::AddEmModel(G4int o <<  99 void G4VMultipleScattering::BuildPhysicsTable(const G4ParticleDefinition& part)
119                                        const G << 
120 {                                                 100 {
121   if(nullptr == ptr) { return; }               << 101   currentCouple = 0;
122   G4VEmFluctuationModel* fm = nullptr;         << 102   if(0 < verboseLevel) {
123   modelManager->AddEmModel(order, ptr, fm, reg << 103     G4cout << "G4VMultipleScattering::BuildPhysicsTable() for "
124   ptr->SetParticleChange(pParticleChange);     << 104            << GetProcessName()
125 }                                              << 105            << " and particle " << part.GetParticleName()
                                                   >> 106            << G4endl;
                                                   >> 107   }
                                                   >> 108 
                                                   >> 109   G4bool cutsWasModified = false;
                                                   >> 110   const G4ProductionCutsTable* theCoupleTable=
                                                   >> 111         G4ProductionCutsTable::GetProductionCutsTable();
                                                   >> 112   size_t numOfCouples = theCoupleTable->GetTableSize();
                                                   >> 113 
                                                   >> 114   for (size_t j=0; j<numOfCouples; j++){
                                                   >> 115     if (theCoupleTable->GetMaterialCutsCouple(j)->IsRecalcNeeded()) {
                                                   >> 116       cutsWasModified = true;
                                                   >> 117       break;
                                                   >> 118     }
                                                   >> 119   }
                                                   >> 120   if( !cutsWasModified ) return;
                                                   >> 121   InitialiseProcess(part);
126                                                   122 
127 //....oooOO0OOooo........oooOO0OOooo........oo << 123   if(latDisplasment) navigator = G4TransportationManager::GetTransportationManager()
                                                   >> 124              ->GetNavigatorForTracking();
                                                   >> 125   modelManager->Clear();
                                                   >> 126   const G4DataVector* theCuts = modelManager->Initialise(&part, 0, 10.0, verboseLevel);
128                                                   127 
129 void G4VMultipleScattering::SetEmModel(G4VMscM << 128   if (buildLambdaTable) {
130 {                                              << 
131   if(nullptr == ptr) { return; }               << 
132   if(!mscModels.empty()) {                     << 
133     for(auto & msc : mscModels) { if(msc == pt << 
134   }                                            << 
135   mscModels.push_back(ptr);                    << 
136 }                                              << 
137                                                   129 
138 //....oooOO0OOooo........oooOO0OOooo........oo << 130     theLambdaTable = new G4PhysicsTable(numOfCouples);
139                                                   131 
140 void                                           << 132     for (size_t i=0; i<numOfCouples; i++) {
141 G4VMultipleScattering::PreparePhysicsTable(con << 
142 {                                              << 
143   G4bool master = emManager->IsMaster();       << 
144   if (nullptr == firstParticle) { firstParticl << 
145                                                   133 
146   emManager->PreparePhysicsTable(&part, this); << 134       // create physics vector and fill it
147   currParticle = nullptr;                      << 135       const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i);
                                                   >> 136       G4PhysicsVector* aVector = PhysicsVector(couple);
                                                   >> 137       modelManager->FillLambdaVector(aVector, couple, false);
148                                                   138 
149   if(firstParticle == &part) {                 << 139       // Insert vector for this material into the table
150     baseMat = emManager->GetTableBuilder()->Ge << 140       theLambdaTable->insert(aVector) ;
151     G4EmTableUtil::PrepareMscProcess(this, par << 
152              stepLimit, facrange,              << 
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     }                                             141     }
164   }                                            << 
165 }                                              << 
166                                                   142 
167 //....oooOO0OOooo........oooOO0OOooo........oo << 143     if(0 < verboseLevel) {
168                                                << 144       G4cout << "Lambda table is built for "
169 void G4VMultipleScattering::BuildPhysicsTable( << 145              << part.GetParticleName()
170 {                                              << 146              << G4endl;
171   G4bool master = emManager->IsMaster();       << 147     }
                                                   >> 148     if(2 < verboseLevel) G4cout << *theLambdaTable << G4endl;
                                                   >> 149     if(5 < verboseLevel) G4cout << theCuts << G4endl;
172                                                   150 
173   if(firstParticle == &part) {                 << 
174     emManager->BuildPhysicsTable(&part);       << 
175   }                                            << 
176   const G4VMultipleScattering* ptr = this;     << 
177   if(!master) {                                << 
178     ptr = static_cast<const G4VMultipleScatter << 
179   }                                               151   }
180                                                   152 
181   G4EmTableUtil::BuildMscProcess(this, ptr, pa << 153   G4String num = part.GetParticleName();
182          numberOfModels, master);              << 154   if (verboseLevel>0 || num == "e-" || num == "mu+" || num == "proton")
                                                   >> 155            PrintInfoDefinition();
                                                   >> 156 
                                                   >> 157   if(0 < verboseLevel) {
                                                   >> 158     G4cout << "G4VMultipleScattering::BuildPhysicsTable() done for "
                                                   >> 159            << GetProcessName()
                                                   >> 160            << " and particle " << part.GetParticleName()
                                                   >> 161            << G4endl;
                                                   >> 162   }
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() const
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 " << G4endl;
397                                 currentSafety) << 244   if (theLambdaTable) {
                                                   >> 245     G4cout << "      Lambda tables from "
                                                   >> 246            << G4BestUnit(MinKinEnergy(),"Energy")
                                                   >> 247      << " to "
                                                   >> 248            << G4BestUnit(MaxKinEnergy(),"Energy")
                                                   >> 249            << " in " << nBins << " bins."
                                                   >> 250            << G4endl;
                                                   >> 251   }
                                                   >> 252   if (1 < verboseLevel) {
                                                   >> 253       G4cout << "LambdaTable address= " << theLambdaTable << G4endl;
                                                   >> 254       if(theLambdaTable) G4cout << (*theLambdaTable) << G4endl;
                                                   >> 255   }
398 }                                                 256 }
399                                                   257 
400 //....oooOO0OOooo........oooOO0OOooo........oo    258 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
401                                                   259 
402 G4double G4VMultipleScattering::GetMeanFreePat << 260 G4PhysicsVector* G4VMultipleScattering::PhysicsVector(const G4MaterialCutsCouple* couple)
403               const G4Track&, G4double, G4Forc << 
404 {                                                 261 {
405   *condition = Forced;                         << 262   G4int nbins = 3;
406   return DBL_MAX;                              << 263   //G4int nbins =  nDEDXBins;
                                                   >> 264   if( couple->IsUsed() ) nbins = nBins;
                                                   >> 265   //  G4double xmax = maxKinEnergy*exp( log(maxKinEnergy/minKinEnergy) / ((G4double)(nbins-1)) );
                                                   >> 266   G4PhysicsVector* v = new G4PhysicsLogVector(minKinEnergy, maxKinEnergy, nbins);
                                                   >> 267   return v;
407 }                                                 268 }
408                                                   269 
409 //....oooOO0OOooo........oooOO0OOooo........oo << 270 G4bool G4VMultipleScattering::StorePhysicsTable(G4ParticleDefinition* part,
410                                                << 271              const G4String& directory,
411 G4bool                                         << 272                    G4bool ascii)
412 G4VMultipleScattering::StorePhysicsTable(const << 
413                                          const << 
414                                          G4boo << 
415 {                                                 273 {
416   G4bool yes = true;                           << 274   G4bool res = true;
417   if(part != firstParticle || !emManager->IsMa << 275   if ( theLambdaTable ) {
                                                   >> 276     const G4String name = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
                                                   >> 277     G4bool yes = theLambdaTable->StorePhysicsTable(name,ascii);
                                                   >> 278     if( !yes ) res = false;
                                                   >> 279   }
418                                                   280 
419   return G4EmTableUtil::StoreMscTable(this, pa << 281   if ( res ) {
420               numberOfModels, verboseLevel,    << 282     G4cout << "Physics table are stored for " << part->GetParticleName()
421                                       ascii);  << 283            << " and process " << GetProcessName()
                                                   >> 284      << " in the directory <" << directory
                                                   >> 285      << "> " << G4endl;
                                                   >> 286   } else {
                                                   >> 287     G4cout << "Fail to store Physics Table for " << part->GetParticleName()
                                                   >> 288            << " and process " << GetProcessName()
                                                   >> 289      << " in the directory <" << directory
                                                   >> 290      << "> " << G4endl;
                                                   >> 291   }
                                                   >> 292   return res;
422 }                                                 293 }
423                                                   294 
424 //....oooOO0OOooo........oooOO0OOooo........oo    295 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
425                                                   296 
426 G4bool                                         << 297 G4bool G4VMultipleScattering::RetrievePhysicsTable(G4ParticleDefinition* part,
427 G4VMultipleScattering::RetrievePhysicsTable(co << 298                   const G4String& directory,
428                                             co << 299                         G4bool ascii)
429                                             G4 << 300 {
430 {                                              << 301   if(0 < verboseLevel) {
431   return true;                                 << 302     G4cout << "G4VMultipleScattering::RetrievePhysicsTable() for "
432 }                                              << 303            << part->GetParticleName() << " and process "
433                                                << 304      << GetProcessName() << G4endl;
434 //....oooOO0OOooo........oooOO0OOooo........oo << 305   }
435                                                << 306   G4bool res = true;
436 void G4VMultipleScattering::ProcessDescription << 307   if(!buildLambdaTable) return res;
437 {                                              << 308 
438   if(nullptr != firstParticle) {               << 309   G4String num = part->GetParticleName();
439     StreamInfo(outFile, *firstParticle, true); << 310   const G4ProductionCutsTable* theCoupleTable=
                                                   >> 311           G4ProductionCutsTable::GetProductionCutsTable();
                                                   >> 312   size_t numOfCouples = theCoupleTable->GetTableSize();
                                                   >> 313 
                                                   >> 314   G4String filename = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
                                                   >> 315   theLambdaTable = new G4PhysicsTable(numOfCouples);
                                                   >> 316   res = theLambdaTable->RetrievePhysicsTable(filename,ascii);
                                                   >> 317   if ( res ) {
                                                   >> 318     if (0 < verboseLevel) {
                                                   >> 319         G4cout << "Lambda table for " << num << " is retrieved from <"
                                                   >> 320                << filename << ">"
                                                   >> 321                << G4endl;
                                                   >> 322     }
                                                   >> 323   } else {
                                                   >> 324     theLambdaTable->clearAndDestroy();
                                                   >> 325     theLambdaTable = 0;
                                                   >> 326     if (0 < verboseLevel) {
                                                   >> 327         G4cout << "Lambda table for " << num << " in file <"
                                                   >> 328                << filename << "> is not exist"
                                                   >> 329                << G4endl;
                                                   >> 330     }
440   }                                               331   }
                                                   >> 332 
                                                   >> 333   if (verboseLevel>0 || num == "e-" || num == "mu+" || num == "proton")
                                                   >> 334            PrintInfoDefinition();
                                                   >> 335   return res;
441 }                                                 336 }
442                                                   337 
443 //....oooOO0OOooo........oooOO0OOooo........oo    338 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
444                                                   339 
445                                                   340