Geant4  10.02.p01
G4VMultipleScattering.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4VMultipleScattering.cc 94120 2015-11-06 09:18:03Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4VMultipleScattering
34 //
35 // Author: Vladimir Ivanchenko on base of Laszlo Urban code
36 //
37 // Creation date: 25.03.2003
38 //
39 // Modifications:
40 //
41 // 13.04.03 Change printout (V.Ivanchenko)
42 // 04-06-03 Fix compilation warnings (V.Ivanchenko)
43 // 16-07-03 Use G4VMscModel interface (V.Ivanchenko)
44 // 03-11-03 Fix initialisation problem in RetrievePhysicsTable (V.Ivanchenko)
45 // 04-11-03 Update PrintInfoDefinition (V.Ivanchenko)
46 // 01-03-04 SampleCosineTheta signature changed
47 // 22-04-04 SampleCosineTheta signature changed back to original
48 // 27-08-04 Add InitialiseForRun method (V.Ivanchneko)
49 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
50 // 11-03-05 Shift verbose level by 1 (V.Ivantchenko)
51 // 15-04-05 optimize internal interface (V.Ivanchenko)
52 // 15-04-05 remove boundary flag (V.Ivanchenko)
53 // 27-10-05 introduce virtual function MscStepLimitation() (V.Ivanchenko)
54 // 12-04-07 Add verbosity at destruction (V.Ivanchenko)
55 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
56 // 11-03-08 Set skin value does not effect step limit type (V.Ivanchenko)
57 // 24-06-09 Removed hidden bin in G4PhysicsVector (V.Ivanchenko)
58 // 04-06-13 Adoptation to MT mode (V.Ivanchenko)
59 //
60 // Class Description:
61 //
62 // It is the generic process of multiple scattering it includes common
63 // part of calculations for all charged particles
64 
65 // -------------------------------------------------------------------
66 //
67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
69 
70 #include "G4VMultipleScattering.hh"
71 #include "G4PhysicalConstants.hh"
72 #include "G4SystemOfUnits.hh"
73 #include "G4LossTableManager.hh"
74 #include "G4MaterialCutsCouple.hh"
75 #include "G4Step.hh"
76 #include "G4ParticleDefinition.hh"
77 #include "G4VEmFluctuationModel.hh"
78 #include "G4UnitsTable.hh"
79 #include "G4ProductionCutsTable.hh"
80 #include "G4Electron.hh"
81 #include "G4GenericIon.hh"
83 #include "G4SafetyHelper.hh"
84 #include "G4ParticleTable.hh"
85 #include "G4ProcessVector.hh"
86 #include "G4ProcessManager.hh"
87 #include <iostream>
88 
89 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
90 
91 static const G4double minSafety = 1.20*CLHEP::nm;
92 static const G4double geomMin = 0.05*CLHEP::nm;
94 
98  numberOfModels(0),
99  firstParticle(nullptr),
100  currParticle(nullptr),
101  stepLimit(fUseSafety),
102  facrange(0.04),
103  latDisplacement(true),
104  isIon(false),
105  fDispBeyondSafety(false),
106  fNewPosition(0.,0.,0.),
107  fNewDirection(0.,0.,1.)
108 {
110  SetVerboseLevel(1);
112  if("ionmsc" == name) { firstParticle = G4GenericIon::GenericIon(); }
113 
115 
117  fIonisation = nullptr;
118 
120  safetyHelper = nullptr;
121  fPositionChanged = false;
122  isActive = false;
123 
124  currentModel = nullptr;
127  emManager->Register(this);
128 }
129 
130 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
131 
133 {
134  /*
135  if(1 < verboseLevel) {
136  G4cout << "G4VMultipleScattering destruct " << GetProcessName()
137  << G4endl;
138  }
139  */
140  delete modelManager;
141  emManager->DeRegister(this);
142 }
143 
144 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
145 
147  const G4Region* region)
148 {
149  G4VEmFluctuationModel* fm = nullptr;
150  modelManager->AddEmModel(order, p, fm, region);
151  if(p) { p->SetParticleChange(pParticleChange); }
152 }
153 
154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
155 
157 {
158  G4int n = mscModels.size();
159  if(index >= n) { for(G4int i=n; i<=index; ++i) { mscModels.push_back(0); } }
160  mscModels[index] = p;
161 }
162 
163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
164 
166 {
167  G4VMscModel* p = nullptr;
168  if(index >= 0 && index < G4int(mscModels.size())) { p = mscModels[index]; }
169  return p;
170 }
171 
172 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
173 
174 G4VEmModel*
176 {
177  return modelManager->GetModel(idx, ver);
178 }
179 
180 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
181 
182 void
184 {
185  G4bool master = true;
186  const G4VMultipleScattering* masterProc =
187  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
188  if(masterProc && masterProc != this) { master = false; }
189 
190  if(!firstParticle) { firstParticle = &part; }
191  if(part.GetParticleType() == "nucleus") {
193  latDisplacement = false;
194  facrange = 0.2;
195  G4String pname = part.GetParticleName();
196  if(pname != "deuteron" && pname != "triton" &&
197  pname != "alpha+" && pname != "helium" &&
198  pname != "alpha" && pname != "He3" &&
199  pname != "hydrogen") {
200 
201  const G4ParticleDefinition* theGenericIon =
203  if(&part == theGenericIon) { isIon = true; }
204 
205  if(theGenericIon && firstParticle != theGenericIon) {
206  G4ProcessManager* pm = theGenericIon->GetProcessManager();
208  size_t n = v->size();
209  for(size_t j=0; j<n; ++j) {
210  if((*v)[j] == this) {
211  firstParticle = theGenericIon;
212  isIon = true;
213  break;
214  }
215  }
216  }
217  }
218  }
219 
220  emManager->PreparePhysicsTable(&part, this, master);
221  currParticle = 0;
222 
223  if(1 < verboseLevel) {
224  G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
225  << GetProcessName()
226  << " and particle " << part.GetParticleName()
227  << " local particle " << firstParticle->GetParticleName()
228  << " isIon= " << isIon
229  << G4endl;
230  }
231 
232  if(firstParticle == &part) {
233 
234  // initialise process
236 
237  // heavy particles and not ions
238  if(!isIon) {
239  if(part.GetPDGMass() > MeV) {
243  } else {
247  }
248  if(latDisplacement) {
250  }
251  }
252  if(master) { SetVerboseLevel(theParameters->Verbose()); }
254 
255  // initialisation of models
257  for(G4int i=0; i<numberOfModels; ++i) {
258  G4VMscModel* msc = static_cast<G4VMscModel*>(modelManager->GetModel(i));
259  msc->SetIonisation(0, firstParticle);
260  msc->SetMasterThread(master);
261  if(0 == i) { currentModel = msc; }
264  msc->SetSkin(theParameters->MscSkin());
265  msc->SetRangeFactor(facrange);
268  G4double emax =
270  msc->SetHighEnergyLimit(emax);
271  }
272 
274  10.0, verboseLevel);
275 
276  if(!safetyHelper) {
278  ->GetSafetyHelper();
280  }
281  }
282 }
283 
284 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
285 
287 {
288  G4String num = part.GetParticleName();
289  if(1 < verboseLevel) {
290  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
291  << GetProcessName()
292  << " and particle " << num
293  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
294  << G4endl;
295  }
296  G4bool master = true;
297  const G4VMultipleScattering* masterProcess =
298  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
299  if(masterProcess && masterProcess != this) { master = false; }
300 
301  if(firstParticle == &part) {
302  /*
303  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
304  << GetProcessName()
305  << " and particle " << num
306  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
307  << " " << this
308  << G4endl;
309  */
311 
312  if(!master) {
313  // initialisation of models
314  G4bool printing = true;
316  /*
317  G4cout << "### G4VMultipleScattering::SlaveBuildPhysicsTable() for "
318  << GetProcessName()
319  << " and particle " << num
320  << " Nmod= " << numberOfModels << " " << this
321  << G4endl;
322  */
323  for(G4int i=0; i<numberOfModels; ++i) {
324  G4VMscModel* msc =
325  static_cast<G4VMscModel*>(GetModelByIndex(i, printing));
326  G4VMscModel* msc0=
327  static_cast<G4VMscModel*>(masterProcess->GetModelByIndex(i,printing));
328  msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false);
329  msc->InitialiseLocal(firstParticle, msc0);
330  }
331  }
332  }
333 
334  // explicitly defined printout by particle name
335  if(1 < verboseLevel ||
336  (0 < verboseLevel && (num == "e-" ||
337  num == "e+" || num == "mu+" ||
338  num == "mu-" || num == "proton"||
339  num == "pi+" || num == "pi-" ||
340  num == "kaon+" || num == "kaon-" ||
341  num == "alpha" || num == "anti_proton" ||
342  num == "GenericIon")))
343  {
344  G4cout << G4endl << GetProcessName()
345  << ": for " << num
346  << " SubType= " << GetProcessSubType()
347  << G4endl;
348  PrintInfo();
350  }
351 
352  if(1 < verboseLevel) {
353  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
354  << GetProcessName()
355  << " and particle " << num
356  << G4endl;
357  }
358 }
359 
360 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
361 
363 {
364  if (0 < verboseLevel) {
365  G4cout << G4endl << GetProcessName()
366  << ": for " << firstParticle->GetParticleName()
367  << " SubType= " << GetProcessSubType()
368  << G4endl;
369  PrintInfo();
371  }
372 }
373 
374 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
375 
377 {
378  G4VEnergyLossProcess* eloss = nullptr;
379  if(track->GetParticleDefinition() != currParticle) {
382  eloss = fIonisation;
383  }
384  /*
385  G4cout << "G4VMultipleScattering::StartTracking Nmod= " << numberOfModels
386  << " " << currParticle->GetParticleName()
387  << " E(MeV)= " << track->GetKineticEnergy()
388  << " Ion= " << eloss << " " << fIonisation << " IsMaster= "
389  << G4LossTableManager::Instance()->IsMaster()
390  << G4endl;
391  */
392  // one model
393  if(1 == numberOfModels) {
394  currentModel->StartTracking(track);
396 
397  // many models
398  } else {
399  for(G4int i=0; i<numberOfModels; ++i) {
400  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i,true));
401  msc->StartTracking(track);
402  if(eloss) { msc->SetIonisation(fIonisation, currParticle); }
403  }
404  }
405 }
406 
407 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
408 
410  const G4Track& track,
411  G4double,
412  G4double currentMinimalStep,
413  G4double&,
414  G4GPILSelection* selection)
415 {
416  // get Step limit proposed by the process
417  *selection = NotCandidateForSelection;
418  physStepLimit = gPathLength = tPathLength = currentMinimalStep;
419 
420  G4double ekin = track.GetKineticEnergy();
421  /*
422  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
423  << " " << currParticle->GetParticleName()
424  << " currMod " << currentModel
425  << G4endl;
426  */
427  // isIon flag is used only to select a model
428  if(isIon) {
429  ekin *= proton_mass_c2/track.GetParticleDefinition()->GetPDGMass();
430  }
431 
432  // select new model
433  if(1 < numberOfModels) {
434  currentModel = static_cast<G4VMscModel*>(
435  SelectModel(ekin,track.GetMaterialCutsCouple()->GetIndex()));
436  }
437  // msc is active is model is active, energy above the limit,
438  // and step size is above the limit;
439  // if it is active msc may limit the step
441  && ekin >= lowestKinEnergy) {
442  isActive = true;
443  tPathLength =
445  if (tPathLength < physStepLimit) {
446  *selection = CandidateForSelection;
447  }
448  } else { isActive = false; }
449 
450  //if(currParticle->GetPDGMass() > GeV)
451  /*
452  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
453  << " " << currParticle->GetParticleName()
454  << " gPathLength= " << gPathLength
455  << " tPathLength= " << tPathLength
456  << " currentMinimalStep= " << currentMinimalStep
457  << " isActive " << isActive << G4endl;
458  */
459  return gPathLength;
460 }
461 
462 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
463 
464 G4double
467 {
468  *condition = Forced;
469  //*condition = NotForced;
470  return DBL_MAX;
471 }
472 
473 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
474 
477 {
482  fPositionChanged = false;
483 
484  G4double geomLength = step.GetStepLength();
485 
486  // very small step - no msc
487  if(!isActive) {
488  tPathLength = geomLength;
489 
490  // sample msc
491  } else {
492  G4double range =
494  track.GetMaterialCutsCouple());
495 
497 
498  // protection against wrong t->g->t conversion
499  /*
500  if(currParticle->GetPDGMass() > 0.9*GeV)
501  G4cout << "G4VMsc::AlongStepDoIt: GeomLength= "
502  << geomLength
503  << " tPathLength= " << tPathLength
504  << " physStepLimit= " << physStepLimit
505  << " dr= " << range - tPathLength
506  << " ekin= " << track.GetKineticEnergy() << G4endl;
507  */
509 
510  // do not sample scattering at the last or at a small step
511  if(tPathLength < range && tPathLength > geomMin) {
512 
515 
516  G4double r2 = displacement.mag2();
517  //G4cout << " R= " << sqrt(r2) << " Rmin= " << sqrt(minDisplacement2)
518  // << " flag= " << fDispBeyondSafety << G4endl;
519  if(r2 > minDisplacement2) {
520 
521  fPositionChanged = true;
522  const G4double sFact = 0.99;
523  G4double postSafety =
524  sFact*(step.GetPreStepPoint()->GetSafety() - geomLength);
525  //G4cout<<" R= "<<sqrt(r2)<<" postSafety= "<<postSafety<<G4endl;
526 
527  // far away from geometry boundary
528  if(postSafety > 0.0 && r2 <= postSafety*postSafety) {
529  fNewPosition += displacement;
530 
531  //near the boundary
532  } else {
533  G4double dispR = std::sqrt(r2);
534  postSafety =
535  sFact*safetyHelper->ComputeSafety(fNewPosition, dispR);
536 
537  // displaced point is definitely within the volume
538  //G4cout<<" R= "<<dispR<<" postSafety= "<<postSafety<<G4endl;
539  if(dispR < postSafety) {
540  fNewPosition += displacement;
541 
542  // optional extra mechanism is applied only if a particle
543  // is stopped by the boundary
544  } else if(fDispBeyondSafety && 0.0 == postSafety) {
545  fNewPosition += displacement;
546  G4double maxshift =
547  std::min(2.0*dispR, geomLength*(physStepLimit/tPathLength - 1.0));
548  G4double dist = 0.0;
549  G4double safety = postSafety + dispR;
551  /*
552  G4cout << "##MSC before Recheck maxshift= " << maxshift
553  << " postsafety= " << postSafety
554  << " Ekin= " << track.GetKineticEnergy()
555  << " " << track.GetDefinition()->GetParticleName()
556  << G4endl;
557  */
558  // check if it is possible to shift to the boundary
559  // and the shift is not large
561  fNewDirection, maxshift, &dist, &safety)
562  && std::abs(dist) < maxshift) {
563  /*
564  G4cout << "##MSC after Recheck dist= " << dist
565  << " postsafety= " << postSafety
566  << " t= " << tPathLength
567  << " g= " << geomLength
568  << " p= " << physStepLimit
569  << G4endl;
570  */
571  // shift is positive
572  if(dist >= 0.0) {
573  tPathLength *= (1.0 + dist/geomLength);
574  fNewPosition += dist*fNewDirection;
575 
576  // shift is negative cannot be larger than geomLength
577  } else {
578  maxshift = std::min(maxshift, geomLength);
579  if(0.0 < maxshift + dist) {
580  G4ThreeVector postpoint = step.GetPostStepPoint()->GetPosition();
582  G4double R2 = (postpoint - point).mag2();
583  G4double newdist = dist;
584  // check not more than 10 extra boundaries
585  for(G4int i=0; i<10; ++i) {
586  dist = 0.0;
588  point, fNewDirection, maxshift, &dist, &safety)
589  && std::abs(newdist + dist) < maxshift) {
590  point += dist*fNewDirection;
591  G4double R2new = (postpoint - point).mag2();
592  //G4cout << "Backward i= " << i << " dist= " << dist
593  // << " R2= " << R2new << G4endl;
594  if(dist >= 0.0 || R2new > R2) { break; }
595  R2 = R2new;
596  fNewPosition = point;
597  newdist += dist;
598  } else {
599  break;
600  }
601  }
602  tPathLength *= (1.0 + newdist/geomLength);
603  // shift on boundary is not possible for negative disp
604  } else {
605  fNewPosition += displacement*(postSafety/dispR - 1.0);
606  }
607  }
608  // shift on boundary is not possible for any disp
609  } else {
610  fNewPosition += displacement*(postSafety/dispR - 1.0);
611  }
612  // reduced displacement
613  } else if(postSafety > geomMin) {
614  fNewPosition += displacement*(postSafety/dispR);
615 
616  // very small postSafety
617  } else {
618  fPositionChanged = false;
619  }
620  }
621  //safetyHelper->ReLocateWithinVolume(fNewPosition);
622  }
623  }
624  }
626  //fParticleChange.ProposePosition(fNewPosition);
627  return &fParticleChange;
628 }
629 
630 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
631 
634 {
636 
637  if(fPositionChanged) {
640  }
641 
642  return &fParticleChange;
643 }
644 
645 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
646 
648  const G4Track& track,
649  G4double previousStepSize,
650  G4double currentMinimalStep,
651  G4double& currentSafety)
652 {
654  G4double x = AlongStepGetPhysicalInteractionLength(track,previousStepSize,
655  currentMinimalStep,
656  currentSafety,
657  &selection);
658  return x;
659 }
660 
661 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
662 
664  const G4Track& track,
665  G4double previousStepSize,
666  G4double currentMinimalStep,
667  G4double& currentSafety)
668 {
669  return GetContinuousStepLimit(track,previousStepSize,currentMinimalStep,
670  currentSafety);
671 }
672 
673 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
674 
677 {
678  *condition = Forced;
679  return DBL_MAX;
680 }
681 
682 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
683 
684 G4bool
686  const G4String& directory,
687  G4bool ascii)
688 {
689  G4bool yes = true;
690  if(part != firstParticle) { return yes; }
691  const G4VMultipleScattering* masterProcess =
692  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
693  if(masterProcess && masterProcess != this) { return yes; }
694 
696  static const G4String ss[4] = {"1","2","3","4"};
697  for(G4int i=0; i<nmod; ++i) {
698  G4VEmModel* msc = modelManager->GetModel(i);
699  yes = true;
700  G4PhysicsTable* table = msc->GetCrossSectionTable();
701  if (table) {
702  G4int j = std::min(i,3);
703  G4String name =
704  GetPhysicsTableFileName(part,directory,"LambdaMod"+ss[j],ascii);
705  yes = table->StorePhysicsTable(name,ascii);
706 
707  if ( yes ) {
708  if ( verboseLevel>0 ) {
709  G4cout << "Physics table are stored for "
710  << part->GetParticleName()
711  << " and process " << GetProcessName()
712  << " with a name <" << name << "> " << G4endl;
713  }
714  } else {
715  G4cout << "Fail to store Physics Table for "
716  << part->GetParticleName()
717  << " and process " << GetProcessName()
718  << " in the directory <" << directory
719  << "> " << G4endl;
720  }
721  }
722  }
723  return yes;
724 }
725 
726 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
727 
728 G4bool
730  const G4String&,
731  G4bool)
732 {
733  return true;
734 }
735 
736 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
737 
739 {
740  for(G4int i=0; i<numberOfModels; ++i) {
741  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i, true));
742  msc->SetIonisation(p, firstParticle);
743  }
744 }
745 
746 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
747 
748 void G4VMultipleScattering::ProcessDescription(std::ostream& outFile) const
749 {
750  outFile << "Multiple scattering process <" << GetProcessName()
751  << ">" << G4endl;
752 }
753 
754 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
755 
G4double condition(const G4ErrorSymMatrix &m)
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii)
const G4ThreeVector * GetMomentumDirection() const
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
static const G4double geomMin
static const double MeV
Definition: G4SIunits.hh:211
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
virtual void Initialize(const G4Track &)
G4SafetyHelper * GetSafetyHelper() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
void ReLocateWithinVolume(const G4ThreeVector &pGlobalPoint)
virtual void StartTracking(G4Track *)
Definition: G4VEmModel.cc:284
void InitialiseHelper()
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
static G4LossTableManager * Instance()
G4EmModelManager * modelManager
G4int verboseLevel
Definition: G4VProcess.hh:368
CLHEP::Hep3Vector G4ThreeVector
G4double GetStepLength() const
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:634
G4MscStepLimitType MscMuHadStepLimitType() const
G4VEmModel * GetModel(G4int, G4bool ver=false)
void DeRegister(G4VEnergyLossProcess *p)
G4double MscGeomFactor() const
G4String name
Definition: TRTMaterials.hh:40
G4double MscMuHadRangeFactor() const
G4double MscThetaLimit() const
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:219
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
virtual void ProcessDescription(std::ostream &outFile) const
G4double GetMeanFreePath(const G4Track &track, G4double, G4ForceCondition *condition)
void SetLateralDisplasmentFlag(G4bool val)
Definition: G4VMscModel.hh:197
G4double PostStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4ForceCondition *condition)
void SetEmModel(G4VMscModel *, G4int index=1)
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:826
void SetStepLimitType(G4MscStepLimitType)
Definition: G4VMscModel.hh:225
G4ProcessManager * GetProcessManager() const
int G4int
Definition: G4Types.hh:78
G4VEnergyLossProcess * fIonisation
const G4String & GetPhysicsTableFileName(const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
Definition: G4VProcess.cc:186
static const G4double minDisplacement2
const G4String & GetParticleName() const
G4bool LatDisplacementBeyondSafety() const
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:725
G4VMultipleScattering(const G4String &name="msc", G4ProcessType type=fElectromagnetic)
G4StepPoint * GetPreStepPoint() const
const G4ParticleDefinition * currParticle
G4ParticleChangeForMSC fParticleChange
virtual G4double ComputeTrueStepLength(G4double geomPathLength)
Definition: G4VMscModel.cc:160
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4GLOB_DLL std::ostream G4cout
G4int Verbose() const
virtual G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &stepLimit)
Definition: G4VMscModel.cc:146
G4bool MuHadLateralDisplacement() const
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=0)
Definition: G4VEmModel.cc:410
const G4ThreeVector & GetPosition() const
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:295
bool G4bool
Definition: G4Types.hh:79
static const double nm
Definition: G4SIunits.hh:111
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
void SetRangeFactor(G4double)
Definition: G4VMscModel.hh:211
void ProposeTrueStepLength(G4double truePathLength)
G4bool LateralDisplacement() const
const G4ParticleDefinition * GetParticleDefinition() const
G4bool IsMaster() const
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
const G4String & GetParticleType() const
void SetCrossSectionTable(G4PhysicsTable *, G4bool isLocal)
Definition: G4VEmModel.cc:418
void Register(G4VEnergyLossProcess *p)
Definition: G4Step.hh:76
const G4int n
G4LossTableManager * emManager
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:753
void DumpModelList(G4int verb)
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
G4double ComputeSafety(const G4ThreeVector &pGlobalPoint, G4double maxRadius=DBL_MAX)
static G4TransportationManager * GetTransportationManager()
static const G4double emax
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:711
static const G4double minSafety
G4int size() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
static const double eV
Definition: G4SIunits.hh:212
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:342
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
G4int NumberOfModels() const
void ProposePosition(const G4ThreeVector &finalPosition)
void BuildPhysicsTable(const G4ParticleDefinition *aParticle)
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
G4VParticleChange * AlongStepDoIt(const G4Track &, const G4Step &)
G4double ContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
const G4ParticleDefinition * firstParticle
G4StepPoint * GetPostStepPoint() const
void AddEmModel(G4int order, G4VEmModel *, const G4Region *region=0)
const G4double x[NPOINTSGL]
static G4EmParameters * Instance()
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double MscRangeFactor() const
G4double GetSafety() const
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
G4MscStepLimitType stepLimit
void SetGeomFactor(G4double)
Definition: G4VMscModel.hh:218
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)
G4MscStepLimitType MscStepLimitType() const
virtual void PrintInfo()=0
double G4double
Definition: G4Types.hh:76
void SetSkin(G4double)
Definition: G4VMscModel.hh:204
void BuildPhysicsTable(const G4ParticleDefinition &)
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4ForceCondition
G4bool RecheckDistanceToCurrentBoundary(const G4ThreeVector &pGlobalPoint, const G4ThreeVector &pDirection, const G4double pCurrentProposedStepLength, G4double *prDistance, G4double *prNewSafety=0) const
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection)
#define DBL_MAX
Definition: templates.hh:83
G4VMscModel * EmModel(G4int index=1) const
G4double MscSkin() const
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426
std::vector< G4VMscModel * > mscModels
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:760
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
G4GPILSelection
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
void PreparePhysicsTable(const G4ParticleDefinition &)
G4ProcessType
virtual G4ThreeVector & SampleScattering(const G4ThreeVector &, G4double safety)
Definition: G4VMscModel.cc:139
void SetIonisation(G4VEnergyLossProcess *)