Geant4  10.01.p03
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 90579 2015-06-04 10:00:26Z 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 
88 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
89 
90 static const G4double minSafety = 1.20*CLHEP::nm;
91 static const G4double geomMin = 0.05*CLHEP::nm;
93 
97  numberOfModels(0),
98  firstParticle(0),
99  currParticle(0),
100  stepLimit(fUseSafety),
101  facrange(0.04),
102  latDisplacement(true),
103  isIon(false),
104  fNewPosition(0.,0.,0.),
105  fNewDirection(0.,0.,1.),
106  fDispBeyondSafety(false)
107 {
109  SetVerboseLevel(1);
111  if("ionmsc" == name) { firstParticle = G4GenericIon::GenericIon(); }
112 
114 
116  fIonisation = 0;
117 
119  safetyHelper = 0;
120  fPositionChanged = false;
121  isActive = false;
122  actStepLimit = false;
123  actFacRange = false;
124  actLatDisp = false;
125 
126  currentModel = 0;
129  emManager->Register(this);
130 }
131 
132 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
133 
135 {
136  /*
137  if(1 < verboseLevel) {
138  G4cout << "G4VMultipleScattering destruct " << GetProcessName()
139  << G4endl;
140  }
141  */
142  delete modelManager;
143  emManager->DeRegister(this);
144 }
145 
146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
147 
149  const G4Region* region)
150 {
151  G4VEmFluctuationModel* fm = 0;
152  modelManager->AddEmModel(order, p, fm, region);
153  if(p) { p->SetParticleChange(pParticleChange); }
154 }
155 
156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
157 
159 {
160  G4int n = mscModels.size();
161  if(index >= n) { for(G4int i=n; i<=index; ++i) { mscModels.push_back(0); } }
162  mscModels[index] = p;
163 }
164 
165 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
166 
168 {
169  G4VMscModel* p = 0;
170  if(index >= 0 && index < G4int(mscModels.size())) { p = mscModels[index]; }
171  return p;
172 }
173 
174 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
175 
176 G4VEmModel*
178 {
179  return modelManager->GetModel(idx, ver);
180 }
181 
182 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
183 
184 void
186 {
187  G4bool master = true;
188  const G4VMultipleScattering* masterProc =
189  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
190  if(masterProc && masterProc != this) { master = false; }
191 
192  if(!firstParticle) { firstParticle = &part; }
193  if(part.GetParticleType() == "nucleus") {
196  SetRangeFactor(0.2);
197  G4String pname = part.GetParticleName();
198  if(pname != "deuteron" && pname != "triton" &&
199  pname != "alpha+" && pname != "helium" &&
200  pname != "alpha" && pname != "He3" &&
201  pname != "hydrogen") {
202 
203  const G4ParticleDefinition* theGenericIon =
205 
206  if(theGenericIon && firstParticle != theGenericIon) {
207  G4ProcessManager* pm = theGenericIon->GetProcessManager();
209  size_t n = v->size();
210  for(size_t j=0; j<n; ++j) {
211  if((*v)[j] == this) {
212  firstParticle = theGenericIon;
213  isIon = true;
214  break;
215  }
216  }
217  }
218  }
219  }
220 
221  emManager->PreparePhysicsTable(&part, this, master);
222  currParticle = 0;
223 
224  if(1 < verboseLevel) {
225  G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
226  << GetProcessName()
227  << " and particle " << part.GetParticleName()
228  << " local particle " << firstParticle->GetParticleName()
229  << " isIon= " << isIon
230  << G4endl;
231  }
232 
233  if(firstParticle == &part) {
234 
235  // initialise process
237  if(part.GetPDGMass() > MeV) {
238  if(!actStepLimit) { stepLimit = fMinimal; }
239  if(!actFacRange) { facrange = 0.2; }
240  if(!actLatDisp) {
242  }
243  } else {
246  if(!actLatDisp) {
248  }
249  }
250  if(latDisplacement) {
252  }
253  if(master) { SetVerboseLevel(theParameters->Verbose()); }
255 
256  // initialisation of models
258  for(G4int i=0; i<numberOfModels; ++i) {
259  G4VMscModel* msc = static_cast<G4VMscModel*>(modelManager->GetModel(i));
260  msc->SetIonisation(0, firstParticle);
261  msc->SetMasterThread(master);
262  if(0 == i) { currentModel = msc; }
265  msc->SetSkin(theParameters->MscSkin());
266  msc->SetRangeFactor(facrange);
269  G4double emax =
271  msc->SetHighEnergyLimit(emax);
272  }
273 
275  10.0, verboseLevel);
276 
277  if(!safetyHelper) {
279  ->GetSafetyHelper();
281  }
282  }
283 }
284 
285 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
286 
288 {
289  G4String num = part.GetParticleName();
290  if(1 < verboseLevel) {
291  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
292  << GetProcessName()
293  << " and particle " << num
294  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
295  << G4endl;
296  }
297  G4bool master = true;
298  const G4VMultipleScattering* masterProcess =
299  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
300  if(masterProcess && masterProcess != this) { master = false; }
301 
302  if(firstParticle == &part) {
303  /*
304  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
305  << GetProcessName()
306  << " and particle " << num
307  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
308  << " " << this
309  << G4endl;
310  */
312 
313  if(!master) {
314  // initialisation of models
315  G4bool printing = true;
317  /*
318  G4cout << "### G4VMultipleScattering::SlaveBuildPhysicsTable() for "
319  << GetProcessName()
320  << " and particle " << num
321  << " Nmod= " << numberOfModels << " " << this
322  << G4endl;
323  */
324  for(G4int i=0; i<numberOfModels; ++i) {
325  G4VMscModel* msc =
326  static_cast<G4VMscModel*>(GetModelByIndex(i, printing));
327  G4VMscModel* msc0=
328  static_cast<G4VMscModel*>(masterProcess->GetModelByIndex(i,printing));
329  msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false);
330  msc->InitialiseLocal(firstParticle, msc0);
331  }
332  }
333  }
334 
335  // explicitly defined printout by particle name
336  if(1 < verboseLevel ||
337  (0 < verboseLevel && (num == "e-" ||
338  num == "e+" || num == "mu+" ||
339  num == "mu-" || num == "proton"||
340  num == "pi+" || num == "pi-" ||
341  num == "kaon+" || num == "kaon-" ||
342  num == "alpha" || num == "anti_proton" ||
343  num == "GenericIon")))
344  {
345  G4cout << G4endl << GetProcessName()
346  << ": for " << num
347  << " SubType= " << GetProcessSubType()
348  << G4endl;
349  PrintInfo();
351  }
352 
353  if(1 < verboseLevel) {
354  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
355  << GetProcessName()
356  << " and particle " << num
357  << G4endl;
358  }
359 }
360 
361 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
362 
364 {
365  if (0 < verboseLevel) {
366  G4cout << G4endl << GetProcessName()
367  << ": for " << firstParticle->GetParticleName()
368  << " SubType= " << GetProcessSubType()
369  << G4endl;
370  PrintInfo();
372  }
373 }
374 
375 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
376 
378 {
379  G4VEnergyLossProcess* eloss = 0;
380  if(track->GetParticleDefinition() != currParticle) {
383  eloss = fIonisation;
384  }
385  /*
386  G4cout << "G4VMultipleScattering::StartTracking Nmod= " << numberOfModels
387  << " " << currParticle->GetParticleName()
388  << " E(MeV)= " << track->GetKineticEnergy()
389  << " Ion= " << eloss << " " << fIonisation << " IsMaster= "
390  << G4LossTableManager::Instance()->IsMaster()
391  << G4endl;
392  */
393  // one model
394  if(1 == numberOfModels) {
395  currentModel->StartTracking(track);
397 
398  // many models
399  } else {
400  for(G4int i=0; i<numberOfModels; ++i) {
401  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i,true));
402  msc->StartTracking(track);
403  if(eloss) { msc->SetIonisation(fIonisation, currParticle); }
404  }
405  }
406 }
407 
408 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
409 
411  const G4Track& track,
412  G4double,
413  G4double currentMinimalStep,
414  G4double&,
415  G4GPILSelection* selection)
416 {
417  // get Step limit proposed by the process
418  *selection = NotCandidateForSelection;
419  physStepLimit = gPathLength = tPathLength = currentMinimalStep;
420 
421  G4double ekin = track.GetKineticEnergy();
422  /*
423  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
424  << " " << currParticle->GetParticleName()
425  << " currMod " << currentModel
426  << G4endl;
427  */
428  // isIon flag is used only to select a model
429  if(isIon) {
430  ekin *= proton_mass_c2/track.GetParticleDefinition()->GetPDGMass();
431  }
432 
433  // select new model
434  if(1 < numberOfModels) {
435  currentModel = static_cast<G4VMscModel*>(
436  SelectModel(ekin,track.GetMaterialCutsCouple()->GetIndex()));
437  }
438  // msc is active is model is active, energy above the limit,
439  // and step size is above the limit;
440  // if it is active msc may limit the step
442  && ekin >= lowestKinEnergy) {
443  isActive = true;
444  tPathLength =
446  if (tPathLength < physStepLimit) {
447  *selection = CandidateForSelection;
448  }
449  } else { isActive = false; }
450 
451  //if(currParticle->GetPDGMass() > GeV)
452  /*
453  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
454  << " " << currParticle->GetParticleName()
455  << " gPathLength= " << gPathLength
456  << " tPathLength= " << tPathLength
457  << " currentMinimalStep= " << currentMinimalStep
458  << " isActive " << isActive << G4endl;
459  */
460  return gPathLength;
461 }
462 
463 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
464 
465 G4double
468 {
469  *condition = Forced;
470  //*condition = NotForced;
471  return DBL_MAX;
472 }
473 
474 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
475 
478 {
483  fPositionChanged = false;
484 
485  G4double geomLength = step.GetStepLength();
486 
487  // very small step - no msc
488  if(!isActive) {
489  tPathLength = geomLength;
490 
491  // sample msc
492  } else {
493  G4double range =
495  track.GetMaterialCutsCouple());
496 
498 
499  // protection against wrong t->g->t conversion
500  /*
501  if(currParticle->GetPDGMass() > 0.9*GeV)
502  G4cout << "G4VMsc::AlongStepDoIt: GeomLength= "
503  << geomLength
504  << " tPathLength= " << tPathLength
505  << " physStepLimit= " << physStepLimit
506  << " dr= " << range - tPathLength
507  << " ekin= " << track.GetKineticEnergy() << G4endl;
508  */
510 
511  // do not sample scattering at the last or at a small step
512  if(tPathLength < range && tPathLength > geomMin) {
513 
516 
517  G4double r2 = displacement.mag2();
518  //G4cout << " R= " << sqrt(r2) << " Rmin= " << sqrt(minDisplacement2)
519  // << " flag= " << fDispBeyondSafety << G4endl;
520  if(r2 > minDisplacement2) {
521 
522  fPositionChanged = true;
523  const G4double sFact = 0.99;
524  G4double postSafety =
525  sFact*(step.GetPreStepPoint()->GetSafety() - geomLength);
526  //G4cout<<" R= "<<sqrt(r2)<<" postSafety= "<<postSafety<<G4endl;
527 
528  // far away from geometry boundary
529  if(postSafety > 0.0 && r2 <= postSafety*postSafety) {
530  fNewPosition += displacement;
531 
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, physStepLimit-tPathLength);
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
561  maxshift,
562  &dist,
563  &safety))
564  {
565  /*
566  G4cout << "##MSC after Recheck dist= " << dist
567  << " postsafety= " << postSafety
568  << " t= " << tPathLength
569  << " g= " << geomLength
570  << " p= " << physStepLimit
571  << G4endl;
572  */
573  G4double tnew = tPathLength*(1.0 + dist/geomLength);
574  if(tnew >= 0.0 && tnew < physStepLimit) {
575  tPathLength = tnew;
576  fNewPosition += dist*fNewDirection;
577  } else {
578  fNewPosition += displacement*(postSafety/dispR - 1.0);
579  }
580  }
581  else
582  // shift on boundary is not possible
583  {
584  fNewPosition += displacement*(postSafety/dispR - 1.0);
585  }
586  // reduced displacement
587  } else if(postSafety > geomMin) {
588  fNewPosition += displacement*(postSafety/dispR);
589 
590  // very small postSafety
591  } else {
592  fPositionChanged = false;
593  }
594  }
595  //safetyHelper->ReLocateWithinVolume(fNewPosition);
596  }
597  }
598  }
600  //fParticleChange.ProposePosition(fNewPosition);
601  return &fParticleChange;
602 }
603 
604 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
605 
608 {
610 
611  if(fPositionChanged) {
614  }
615 
616  return &fParticleChange;
617 }
618 
619 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
620 
622  const G4Track& track,
623  G4double previousStepSize,
624  G4double currentMinimalStep,
625  G4double& currentSafety)
626 {
628  G4double x = AlongStepGetPhysicalInteractionLength(track,previousStepSize,
629  currentMinimalStep,
630  currentSafety,
631  &selection);
632  return x;
633 }
634 
635 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
636 
638  const G4Track& track,
639  G4double previousStepSize,
640  G4double currentMinimalStep,
641  G4double& currentSafety)
642 {
643  return GetContinuousStepLimit(track,previousStepSize,currentMinimalStep,
644  currentSafety);
645 }
646 
647 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
648 
651 {
652  *condition = Forced;
653  return DBL_MAX;
654 }
655 
656 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
657 
658 G4bool
660  const G4String& directory,
661  G4bool ascii)
662 {
663  G4bool yes = true;
664  if(part != firstParticle) { return yes; }
665  const G4VMultipleScattering* masterProcess =
666  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
667  if(masterProcess && masterProcess != this) { return yes; }
668 
670  static const G4String ss[4] = {"1","2","3","4"};
671  for(G4int i=0; i<nmod; ++i) {
672  G4VEmModel* msc = modelManager->GetModel(i);
673  yes = true;
674  G4PhysicsTable* table = msc->GetCrossSectionTable();
675  if (table) {
676  G4int j = std::min(i,3);
677  G4String name =
678  GetPhysicsTableFileName(part,directory,"LambdaMod"+ss[j],ascii);
679  yes = table->StorePhysicsTable(name,ascii);
680 
681  if ( yes ) {
682  if ( verboseLevel>0 ) {
683  G4cout << "Physics table are stored for "
684  << part->GetParticleName()
685  << " and process " << GetProcessName()
686  << " with a name <" << name << "> " << G4endl;
687  }
688  } else {
689  G4cout << "Fail to store Physics Table for "
690  << part->GetParticleName()
691  << " and process " << GetProcessName()
692  << " in the directory <" << directory
693  << "> " << G4endl;
694  }
695  }
696  }
697  return yes;
698 }
699 
700 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
701 
702 G4bool
704  const G4String&,
705  G4bool)
706 {
707  return true;
708 }
709 
710 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
711 
713 {
714  for(G4int i=0; i<numberOfModels; ++i) {
715  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i, true));
716  msc->SetIonisation(p, firstParticle);
717  }
718 }
719 
720 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
721 
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:193
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:616
void SetLateralDisplasmentFlag(G4bool val)
G4VEmModel * GetModel(G4int, G4bool ver=false)
void DeRegister(G4VEnergyLossProcess *p)
G4double MscGeomFactor() const
G4String name
Definition: TRTMaterials.hh:40
G4double MscThetaLimit() const
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:219
const G4MaterialCutsCouple * GetMaterialCutsCouple() 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:808
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:707
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:402
const G4ThreeVector & GetPosition() const
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:289
bool G4bool
Definition: G4Types.hh:79
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:410
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:735
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()
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:693
static const G4double minSafety
G4int size() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
static const double eV
Definition: G4SIunits.hh:194
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:336
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)
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
void SetRangeFactor(G4double val)
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:742
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
G4GPILSelection
void SetStepLimitType(G4MscStepLimitType val)
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 *)