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 7.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.30 2004/12/01 18:01:01 vnivanch Exp $
                                                   >>  24 // GEANT4 tag $Name: geant4-07-00-cand-03 $
 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 // 01-03-04 SampleCosineTheta signature changed
 44 // 22-04-04 SampleCosineTheta signature change     45 // 22-04-04 SampleCosineTheta signature changed back to original
 45 // 27-08-04 Add InitialiseForRun method (V.Iva     46 // 27-08-04 Add InitialiseForRun method (V.Ivanchneko)
 46 // 08-11-04 Migration to new interface of Stor     47 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
 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 //                                                 48 //
                                                   >>  49 // Class Description:
                                                   >>  50 //
                                                   >>  51 // It is the generic process of multiple scattering it includes common
                                                   >>  52 // part of calculations for all charged particles
 57                                                    53 
 58 // -------------------------------------------     54 // -------------------------------------------------------------------
 59 //                                                 55 //
 60 //....oooOO0OOooo........oooOO0OOooo........oo     56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 61 //....oooOO0OOooo........oooOO0OOooo........oo     57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 62                                                    58 
 63 #include "G4VMultipleScattering.hh"                59 #include "G4VMultipleScattering.hh"
 64 #include "G4PhysicalConstants.hh"              << 
 65 #include "G4SystemOfUnits.hh"                  << 
 66 #include "G4LossTableManager.hh"                   60 #include "G4LossTableManager.hh"
 67 #include "G4MaterialCutsCouple.hh"             << 
 68 #include "G4Step.hh"                               61 #include "G4Step.hh"
 69 #include "G4ParticleDefinition.hh"                 62 #include "G4ParticleDefinition.hh"
 70 #include "G4VEmFluctuationModel.hh"                63 #include "G4VEmFluctuationModel.hh"
                                                   >>  64 #include "G4DataVector.hh"
                                                   >>  65 #include "G4PhysicsTable.hh"
                                                   >>  66 #include "G4PhysicsVector.hh"
                                                   >>  67 #include "G4PhysicsLogVector.hh"
 71 #include "G4UnitsTable.hh"                         68 #include "G4UnitsTable.hh"
 72 #include "G4ProductionCutsTable.hh"                69 #include "G4ProductionCutsTable.hh"
 73 #include "G4Electron.hh"                       <<  70 #include "G4Region.hh"
 74 #include "G4GenericIon.hh"                     <<  71 #include "G4RegionStore.hh"
                                                   >>  72 #include "G4Navigator.hh"
 75 #include "G4TransportationManager.hh"              73 #include "G4TransportationManager.hh"
 76 #include "G4SafetyHelper.hh"                   <<  74 #include "G4PhysicsTableHelper.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                                                    75 
 84 //....oooOO0OOooo........oooOO0OOooo........oo     76 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
 85                                                    77 
 86 G4VMultipleScattering::G4VMultipleScattering(c <<  78 G4VMultipleScattering::G4VMultipleScattering(const G4String& name, G4ProcessType type):
 87   : G4VContinuousDiscreteProcess("msc", fElect <<  79                  G4VContinuousDiscreteProcess(name, type),
 88   fNewPosition(0.,0.,0.),                      <<  80   navigator(0),
 89   fNewDirection(0.,0.,1.)                      <<  81   theLambdaTable(0),
 90 {                                              <<  82   firstParticle(0),
 91   theParameters = G4EmParameters::Instance();  <<  83   currentParticle(0),
 92   SetVerboseLevel(1);                          <<  84   currentCouple(0),
 93   SetProcessSubType(fMultipleScattering);      <<  85   nBins(120),
 94                                                <<  86   boundary(false),
 95   lowestKinEnergy = 10*CLHEP::eV;              <<  87   latDisplasment(true),
 96                                                <<  88   buildLambdaTable(true)
 97   geomMin = 0.05*CLHEP::nm;                    <<  89 {
 98   minDisplacement2 = geomMin*geomMin;          <<  90   minKinEnergy = 100.0*eV;
 99                                                <<  91   maxKinEnergy = 100.0*TeV;
100   pParticleChange = &fParticleChange;          <<  92   //SetVerboseLevel(0);
101                                                << 
102   modelManager = new G4EmModelManager();           93   modelManager = new G4EmModelManager();
103   emManager = G4LossTableManager::Instance();  <<  94   (G4LossTableManager::Instance())->Register(this);
104   mscModels.reserve(2);                        <<  95 
105   emManager->Register(this);                   << 
106 }                                                  96 }
107                                                    97 
108 //....oooOO0OOooo........oooOO0OOooo........oo     98 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
109                                                    99 
110 G4VMultipleScattering::~G4VMultipleScattering(    100 G4VMultipleScattering::~G4VMultipleScattering()
111 {                                                 101 {
112   delete modelManager;                            102   delete modelManager;
113   emManager->DeRegister(this);                 << 103   if (theLambdaTable) theLambdaTable->clearAndDestroy();
                                                   >> 104   (G4LossTableManager::Instance())->DeRegister(this);
114 }                                                 105 }
115                                                   106 
116 //....oooOO0OOooo........oooOO0OOooo........oo    107 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
117                                                   108 
118 void G4VMultipleScattering::AddEmModel(G4int o << 109 void G4VMultipleScattering::BuildPhysicsTable(const G4ParticleDefinition& part)
119                                        const G << 
120 {                                                 110 {
121   if(nullptr == ptr) { return; }               << 111   if(0 < verboseLevel) {
122   G4VEmFluctuationModel* fm = nullptr;         << 112 //    G4cout << "========================================================" << G4endl;
123   modelManager->AddEmModel(order, ptr, fm, reg << 113     G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
124   ptr->SetParticleChange(pParticleChange);     << 114            << GetProcessName()
125 }                                              << 115            << " and particle " << part.GetParticleName()
                                                   >> 116            << G4endl;
                                                   >> 117   }
126                                                   118 
127 //....oooOO0OOooo........oooOO0OOooo........oo << 119   if (buildLambdaTable && firstParticle == &part) {
128                                                   120 
129 void G4VMultipleScattering::SetEmModel(G4VMscM << 121     const G4ProductionCutsTable* theCoupleTable=
130 {                                              << 122           G4ProductionCutsTable::GetProductionCutsTable();
131   if(nullptr == ptr) { return; }               << 123     size_t numOfCouples = theCoupleTable->GetTableSize();
132   if(!mscModels.empty()) {                     << 
133     for(auto & msc : mscModels) { if(msc == pt << 
134   }                                            << 
135   mscModels.push_back(ptr);                    << 
136 }                                              << 
137                                                   124 
138 //....oooOO0OOooo........oooOO0OOooo........oo << 125     for (size_t i=0; i<numOfCouples; i++) {
139                                                   126 
140 void                                           << 127       if (theLambdaTable->GetFlag(i)) {
141 G4VMultipleScattering::PreparePhysicsTable(con << 
142 {                                              << 
143   G4bool master = emManager->IsMaster();       << 
144   if (nullptr == firstParticle) { firstParticl << 
145                                                   128 
146   emManager->PreparePhysicsTable(&part, this); << 129         // create physics vector and fill it
147   currParticle = nullptr;                      << 130         const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(i);
                                                   >> 131         G4PhysicsVector* aVector = PhysicsVector(couple);
                                                   >> 132         modelManager->FillLambdaVector(aVector, couple, false);
                                                   >> 133         G4PhysicsTableHelper::SetPhysicsVector(theLambdaTable, i, aVector);
                                                   >> 134       }
                                                   >> 135     }
148                                                   136 
149   if(firstParticle == &part) {                 << 137     G4String num = part.GetParticleName();
150     baseMat = emManager->GetTableBuilder()->Ge << 138     if(0 < verboseLevel) {
151     G4EmTableUtil::PrepareMscProcess(this, par << 139       G4cout << "Lambda table is built for "
152              stepLimit, facrange,              << 140              << num
153              latDisplacement, master,          << 141              << G4endl;
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     }                                             142     }
                                                   >> 143     if((verboseLevel>=0) && ( num == "e-" || num == "mu+" || num == "proton" || num == "pi-")) 
                                                   >> 144       PrintInfoDefinition();
                                                   >> 145     if(2 < verboseLevel) G4cout << *theLambdaTable << G4endl;
                                                   >> 146   }
                                                   >> 147 
                                                   >> 148   if(0 < verboseLevel) {
                                                   >> 149     G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
                                                   >> 150            << GetProcessName()
                                                   >> 151            << " and particle " << part.GetParticleName()
                                                   >> 152            << G4endl;
164   }                                               153   }
165 }                                                 154 }
166                                                   155 
167 //....oooOO0OOooo........oooOO0OOooo........oo    156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
168                                                   157 
169 void G4VMultipleScattering::BuildPhysicsTable( << 158 void G4VMultipleScattering::PreparePhysicsTable(const G4ParticleDefinition& part)
170 {                                                 159 {
171   G4bool master = emManager->IsMaster();       << 160   if (!firstParticle) {
172                                                << 161     currentCouple = 0;
173   if(firstParticle == &part) {                 << 162     firstParticle = &part;
174     emManager->BuildPhysicsTable(&part);       << 163     currentParticle = &part;
175   }                                               164   }
176   const G4VMultipleScattering* ptr = this;     << 165 
177   if(!master) {                                << 166   if(0 < verboseLevel) {
178     ptr = static_cast<const G4VMultipleScatter << 167 //    G4cout << "========================================================" << G4endl;
                                                   >> 168     G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
                                                   >> 169            << GetProcessName()
                                                   >> 170            << " and particle " << part.GetParticleName()
                                                   >> 171            << " local particle " << firstParticle->GetParticleName()
                                                   >> 172            << G4endl;
179   }                                               173   }
180                                                   174 
181   G4EmTableUtil::BuildMscProcess(this, ptr, pa << 175   if(firstParticle == &part) {
182          numberOfModels, master);              << 176 
                                                   >> 177     InitialiseProcess(firstParticle);
                                                   >> 178     if(buildLambdaTable)
                                                   >> 179       theLambdaTable = G4PhysicsTableHelper::PreparePhysicsTable(theLambdaTable);
                                                   >> 180     if(latDisplasment && !navigator)
                                                   >> 181       navigator = G4TransportationManager::GetTransportationManager()
                                                   >> 182         ->GetNavigatorForTracking();
                                                   >> 183     const G4DataVector* theCuts = modelManager->Initialise(firstParticle, 0, 10.0, verboseLevel);
                                                   >> 184 
                                                   >> 185     if(2 < verboseLevel) G4cout << theCuts << G4endl;
                                                   >> 186 
                                                   >> 187   }
183 }                                                 188 }
184                                                   189 
185 //....oooOO0OOooo........oooOO0OOooo........oo    190 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
186                                                   191 
187 void G4VMultipleScattering::StreamInfo(std::os << 192 void G4VMultipleScattering::AddEmModel(G4int order, G4VEmModel* p,
188                   const G4ParticleDefinition&  << 193                                  const G4Region* region)
189 {                                                 194 {
190   G4String indent = (rst ? "  " : "");         << 195   G4VEmFluctuationModel* fm = 0;
191   outFile << G4endl << indent << GetProcessNam << 196   modelManager->AddEmModel(order, p, fm, region);
192   if (!rst) outFile << " for " << part.GetPart << 
193   outFile  << "  SubType= " << GetProcessSubTy << 
194   modelManager->DumpModelList(outFile, verbose << 
195 }                                                 197 }
196                                                   198 
197 //....oooOO0OOooo........oooOO0OOooo........oo    199 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
198                                                   200 
199 void G4VMultipleScattering::StartTracking(G4Tr << 201 G4VParticleChange* G4VMultipleScattering::PostStepDoIt(const G4Track& track,
                                                   >> 202                                                        const G4Step& step)
200 {                                                 203 {
201   G4VEnergyLossProcess* eloss = nullptr;       << 204   fParticleChange.Initialize(track);
202   if(track->GetParticleDefinition() != currPar << 205   G4double kineticEnergy = track.GetKineticEnergy();
203     currParticle = track->GetParticleDefinitio << 206   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                                                   207 
216 //....oooOO0OOooo........oooOO0OOooo........oo << 208   if (kineticEnergy > 0.0) {
                                                   >> 209     G4double cth  = currentModel->SampleCosineTheta(truestep,kineticEnergy);
                                                   >> 210     G4double sth  = std::sqrt(1.-cth*cth);
                                                   >> 211     G4double phi  = twopi*G4UniformRand();
                                                   >> 212     G4double dirx = sth*std::cos(phi);
                                                   >> 213     G4double diry = sth*std::sin(phi);
217                                                   214 
218 G4double G4VMultipleScattering::AlongStepGetPh << 215     G4ThreeVector oldDirection = track.GetMomentumDirection();
219                              const G4Track& tr << 216     G4ThreeVector newDirection(dirx,diry,cth);
220                              G4double,         << 217     newDirection.rotateUz(oldDirection);
221                              G4double currentM << 218     fParticleChange.ProposeMomentumDirection(newDirection);
222                              G4double&,        << 
223                              G4GPILSelection*  << 
224 {                                              << 
225   // get Step limit proposed by the process    << 
226   *selection = NotCandidateForSelection;       << 
227   physStepLimit = gPathLength = tPathLength =  << 
228                                                   219 
229   G4double ekin = track.GetKineticEnergy();    << 
230   /*                                              220   /*
231   G4cout << "MSC::AlongStepGPIL: Ekin= " << ek << 221   if(0 < verboseLevel) {
232          << "  " << currParticle->GetParticleN << 222     const G4ParticleDefinition* pd = dynParticle->GetDefinition();
233          << " currMod " << currentModel        << 223     G4cout << "G4VMultipleScattering::PostStepDoIt: Sample secondary; E= " << finalT/MeV
234          << G4endl;                            << 224            << " MeV; model= (" << currentModel->LowEnergyLimit(pd)
235   */                                           << 225            << ", " <<  currentModel->HighEnergyLimit(pd) << ")"
236   // isIon flag is used only to select a model << 226            << 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   }                                               227   }
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   */                                              228   */
273   return gPathLength;                          << 
274 }                                              << 
275                                                   229 
276 //....oooOO0OOooo........oooOO0OOooo........oo << 230     if (latDisplasment) {
277                                                   231 
278 G4double                                       << 232       G4double safety = step.GetPostStepPoint()->GetSafety();
279 G4VMultipleScattering::PostStepGetPhysicalInte << 233       if ( safety > 0.0) {
280               const G4Track&, G4double, G4Forc << 234         G4double r = currentModel->SampleDisplacement();
281 {                                              << 235         if (r > safety) r = safety;
282   *condition = NotForced;                      << 
283   return DBL_MAX;                              << 
284 }                                              << 
285                                                   236 
286 //....oooOO0OOooo........oooOO0OOooo........oo << 237         // sample direction of lateral displacement
                                                   >> 238         G4double phi  = twopi*G4UniformRand();
                                                   >> 239         G4double dirx = std::cos(phi);
                                                   >> 240         G4double diry = std::sin(phi);
287                                                   241 
288 G4VParticleChange*                             << 242         G4ThreeVector latDirection(dirx,diry,0.0);
289 G4VMultipleScattering::AlongStepDoIt(const G4T << 243         latDirection.rotateUz(oldDirection);
290 {                                              << 
291   fParticleChange.InitialiseMSC(track, step);  << 
292   fNewPosition = fParticleChange.GetProposedPo << 
293   fPositionChanged = false;                    << 
294                                                   244 
295   G4double geomLength = step.GetStepLength();  << 245         // compute new endpoint of the Step
                                                   >> 246         G4ThreeVector newPosition = (step.GetPostStepPoint())->GetPosition()
                                                   >> 247               + r*latDirection;
296                                                   248 
297   // very small step - no msc                  << 249         navigator->LocateGlobalPointWithinVolume(newPosition);
298   if(!isActive) {                              << 
299     tPathLength = geomLength;                  << 
300                                                   250 
301     // sample msc                              << 251         fParticleChange.ProposePosition(newPosition);
302   } else {                                     << 
303     G4double range =                           << 
304       currentModel->GetRange(currParticle,trac << 
305                              track.GetMaterial << 
306                                                << 
307     tPathLength = currentModel->ComputeTrueSte << 
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       }                                           252       }
366     }                                             253     }
367   }                                               254   }
368   fParticleChange.ProposeTrueStepLength(tPathL << 255 
369   return &fParticleChange;                        256   return &fParticleChange;
370 }                                                 257 }
371                                                   258 
372 //....oooOO0OOooo........oooOO0OOooo........oo    259 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
373                                                   260 
374 G4double G4VMultipleScattering::GetContinuousS << 261 void G4VMultipleScattering::PrintInfoDefinition()
375                                        const G << 262 {
376                                        G4doubl << 263   G4cout << G4endl << GetProcessName() << ":  Model variant of multiple scattering "
377                                        G4doubl << 264          << "for " << firstParticle->GetParticleName()
378                                        G4doubl << 265          << G4endl;
379 {                                              << 266   if (theLambdaTable) {
380   G4GPILSelection selection = NotCandidateForS << 267     G4cout << "      Lambda tables from "
381   G4double x = AlongStepGetPhysicalInteraction << 268            << G4BestUnit(MinKinEnergy(),"Energy")
382                                                << 269      << " to "
383                                                << 270            << G4BestUnit(MaxKinEnergy(),"Energy")
384                                                << 271            << " in " << nBins << " bins."
385   return x;                                    << 272            << G4endl;
                                                   >> 273   }
                                                   >> 274   if (1 < verboseLevel) {
                                                   >> 275       G4cout << "LambdaTable address= " << theLambdaTable << G4endl;
                                                   >> 276       if(theLambdaTable) G4cout << (*theLambdaTable) << G4endl;
                                                   >> 277   }
386 }                                                 278 }
387                                                   279 
388 //....oooOO0OOooo........oooOO0OOooo........oo    280 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
389                                                   281 
390 G4double G4VMultipleScattering::ContinuousStep << 282 G4PhysicsVector* G4VMultipleScattering::PhysicsVector(const G4MaterialCutsCouple* couple)
391                                        const G << 
392                                        G4doubl << 
393                                        G4doubl << 
394                                        G4doubl << 
395 {                                                 283 {
396   return GetContinuousStepLimit(track,previous << 284   G4int nbins = 3;
397                                 currentSafety) << 285   //G4int nbins =  nDEDXBins;
                                                   >> 286   if( couple->IsUsed() ) nbins = nBins;
                                                   >> 287   //  G4double xmax = maxKinEnergy*std::exp( std::log(maxKinEnergy/minKinEnergy) / ((G4double)(nbins-1)) );
                                                   >> 288   G4PhysicsVector* v = new G4PhysicsLogVector(minKinEnergy, maxKinEnergy, nbins);
                                                   >> 289   return v;
398 }                                                 290 }
399                                                   291 
400 //....oooOO0OOooo........oooOO0OOooo........oo << 292 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
401                                                   293 
402 G4double G4VMultipleScattering::GetMeanFreePat << 294 G4bool G4VMultipleScattering::StorePhysicsTable(const G4ParticleDefinition* part,
403               const G4Track&, G4double, G4Forc << 295                         const G4String& directory,
                                                   >> 296                               G4bool ascii)
404 {                                                 297 {
405   *condition = Forced;                         << 298   G4bool yes = true;
406   return DBL_MAX;                              << 299   if ( theLambdaTable && part == firstParticle) {
                                                   >> 300     const G4String name = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
                                                   >> 301     G4bool yes = theLambdaTable->StorePhysicsTable(name,ascii);
                                                   >> 302 
                                                   >> 303     if ( yes ) {
                                                   >> 304       G4cout << "Physics table are stored for " << part->GetParticleName()
                                                   >> 305              << " and process " << GetProcessName()
                                                   >> 306        << " in the directory <" << directory
                                                   >> 307        << "> " << G4endl;
                                                   >> 308     } else {
                                                   >> 309       G4cout << "Fail to store Physics Table for " << part->GetParticleName()
                                                   >> 310              << " and process " << GetProcessName()
                                                   >> 311        << " in the directory <" << directory
                                                   >> 312        << "> " << G4endl;
                                                   >> 313     }
                                                   >> 314   }
                                                   >> 315   return yes;
407 }                                                 316 }
408                                                   317 
409 //....oooOO0OOooo........oooOO0OOooo........oo    318 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
410                                                   319 
411 G4bool                                         << 320 G4bool G4VMultipleScattering::RetrievePhysicsTable(const G4ParticleDefinition* part,
412 G4VMultipleScattering::StorePhysicsTable(const << 321                      const G4String& directory,
413                                          const << 322                                    G4bool ascii)
414                                          G4boo << 323 {
415 {                                              << 324   if(0 < verboseLevel) {
                                                   >> 325 //    G4cout << "========================================================" << G4endl;
                                                   >> 326     G4cout << "G4VMultipleScattering::RetrievePhysicsTable() for "
                                                   >> 327            << part->GetParticleName() << " and process "
                                                   >> 328      << GetProcessName() << G4endl;
                                                   >> 329   }
416   G4bool yes = true;                              330   G4bool yes = true;
417   if(part != firstParticle || !emManager->IsMa << 
418                                                   331 
419   return G4EmTableUtil::StoreMscTable(this, pa << 332   if(!buildLambdaTable || firstParticle != part) return yes;
420               numberOfModels, verboseLevel,    << 
421                                       ascii);  << 
422 }                                              << 
423                                                   333 
424 //....oooOO0OOooo........oooOO0OOooo........oo << 334   const G4String particleName = part->GetParticleName();
425                                                   335 
426 G4bool                                         << 336   G4String filename = GetPhysicsTableFileName(part,directory,"Lambda",ascii);
427 G4VMultipleScattering::RetrievePhysicsTable(co << 337   yes = G4PhysicsTableHelper::RetrievePhysicsTable(theLambdaTable,filename,ascii);
428                                             co << 338   if ( yes ) {
429                                             G4 << 339     if (-1 < verboseLevel) {
430 {                                              << 340         G4cout << "Lambda table for " << part->GetParticleName() << " is retrieved from <"
431   return true;                                 << 341                << filename << ">"
                                                   >> 342                << G4endl;
                                                   >> 343     }
                                                   >> 344   } else {
                                                   >> 345     if (-1 < verboseLevel) {
                                                   >> 346         G4cout << "Lambda table for " << part->GetParticleName() << " in file <"
                                                   >> 347                << filename << "> is not exist"
                                                   >> 348                << G4endl;
                                                   >> 349     }
                                                   >> 350   }
                                                   >> 351   return yes;
432 }                                                 352 }
433                                                   353 
434 //....oooOO0OOooo........oooOO0OOooo........oo    354 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
435                                                   355 
436 void G4VMultipleScattering::ProcessDescription << 356 const G4PhysicsTable* G4VMultipleScattering::LambdaTable() const
437 {                                                 357 {
438   if(nullptr != firstParticle) {               << 358   return theLambdaTable;
439     StreamInfo(outFile, *firstParticle, true); << 
440   }                                            << 
441 }                                                 359 }
442                                                   360 
443 //....oooOO0OOooo........oooOO0OOooo........oo    361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
444                                                   362 
445                                                   363