Geant4_10
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 78171 2013-12-04 13:22:18Z gunter $
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 
93  numberOfModels(0),
94  firstParticle(0),
95  currParticle(0),
96  stepLimit(fUseSafety),
97  skin(1.0),
98  facrange(0.04),
99  facgeom(2.5),
100  latDisplasment(true),
101  isIon(false)
102 {
103  SetVerboseLevel(1);
105  if("ionmsc" == name) { firstParticle = G4GenericIon::GenericIon(); }
106 
107  geomMin = 0.01*CLHEP::nm;
108  lowestKinEnergy = 10*eV;
109 
110  // default limit on polar angle
111  polarAngleLimit = 0.0;
112 
113  physStepLimit = gPathLength = tPathLength = 0.0;
114  fIonisation = 0;
115 
117  safetyHelper = 0;
118  fPositionChanged = false;
119  isActive = false;
120 
121  modelManager = new G4EmModelManager();
122  emManager = G4LossTableManager::Instance();
123  emManager->Register(this);
124 
125  warn = 0;
126 }
127 
128 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
129 
131 {
132  /*
133  if(1 < verboseLevel) {
134  G4cout << "G4VMultipleScattering destruct " << GetProcessName()
135  << G4endl;
136  }
137  */
138  delete modelManager;
139  emManager->DeRegister(this);
140 }
141 
142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
143 
145  const G4Region* region)
146 {
148  modelManager->AddEmModel(order, p, fm, region);
149  if(p) { p->SetParticleChange(pParticleChange); }
150 }
151 
152 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
153 
155 {
156  G4int n = mscModels.size();
157  if(index >= n) { for(G4int i=n; i<=index; ++i) { mscModels.push_back(0); } }
158  mscModels[index] = p;
159 }
160 
161 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
162 
164 {
165  G4VMscModel* p = 0;
166  if(index >= 0 && index < G4int(mscModels.size())) { p = mscModels[index]; }
167  return p;
168 }
169 
170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
171 
172 G4VEmModel*
174 {
175  return modelManager->GetModel(idx, ver);
176 }
177 
178 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
179 
180 void
182 {
183  G4bool master = true;
184  if(GetMasterProcess() != this) { master = false; }
185 
186  if(!firstParticle) { firstParticle = &part; }
187  if(part.GetParticleType() == "nucleus") {
190  SetRangeFactor(0.2);
191  G4String pname = part.GetParticleName();
192  if(pname != "deuteron" && pname != "triton" &&
193  pname != "alpha+" && pname != "helium" &&
194  pname != "alpha" && pname != "He3" &&
195  pname != "hydrogen") {
196 
197  const G4ParticleDefinition* theGenericIon =
199 
200  if(theGenericIon && firstParticle != theGenericIon) {
201  G4ProcessManager* pm = theGenericIon->GetProcessManager();
203  size_t n = v->size();
204  for(size_t j=0; j<n; ++j) {
205  if((*v)[j] == this) {
206  firstParticle = theGenericIon;
207  isIon = true;
208  break;
209  }
210  }
211  }
212  }
213  }
214 
215  emManager->PreparePhysicsTable(&part, this, master);
216  currParticle = 0;
217 
218  if(1 < verboseLevel) {
219  G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
220  << GetProcessName()
221  << " and particle " << part.GetParticleName()
222  << " local particle " << firstParticle->GetParticleName()
223  << " isIon= " << isIon
224  << G4endl;
225  }
226 
227  if(firstParticle == &part) {
228 
229  InitialiseProcess(firstParticle);
230 
231  // initialisation of models
232  numberOfModels = modelManager->NumberOfModels();
233  for(G4int i=0; i<numberOfModels; ++i) {
234  G4VMscModel* msc = static_cast<G4VMscModel*>(modelManager->GetModel(i));
235  msc->SetIonisation(0, firstParticle);
236  msc->SetMasterThread(master);
237  if(0 == i) { currentModel = msc; }
238  if(isIon) {
240  msc->SetLateralDisplasmentFlag(false);
241  msc->SetRangeFactor(0.2);
242  } else {
245  msc->SetSkin(Skin());
246  msc->SetRangeFactor(RangeFactor());
247  msc->SetGeomFactor(GeomFactor());
248  }
249  msc->SetPolarAngleLimit(polarAngleLimit);
250  G4double emax =
251  std::min(msc->HighEnergyLimit(),emManager->MaxKinEnergy());
252  msc->SetHighEnergyLimit(emax);
253  }
254 
255  modelManager->Initialise(firstParticle, G4Electron::Electron(),
256  10.0, verboseLevel);
257 
258  if(!safetyHelper) {
260  ->GetSafetyHelper();
261  safetyHelper->InitialiseHelper();
262  }
263  }
264 }
265 
266 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
267 
269 {
270  G4String num = part.GetParticleName();
271  if(1 < verboseLevel) {
272  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
273  << GetProcessName()
274  << " and particle " << num
275  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
276  << G4endl;
277  }
278  G4bool master = true;
279  const G4VMultipleScattering* masterProcess =
280  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
281  if(masterProcess != this) { master = false; }
282 
283  if(firstParticle == &part) {
284  /*
285  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
286  << GetProcessName()
287  << " and particle " << num
288  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
289  << " " << this
290  << G4endl;
291  */
292  emManager->BuildPhysicsTable(firstParticle);
293 
294  if(!master) {
295  // initialisation of models
296  G4bool printing = true;
297  numberOfModels = modelManager->NumberOfModels();
298  /*
299  G4cout << "### G4VMultipleScattering::SlaveBuildPhysicsTable() for "
300  << GetProcessName()
301  << " and particle " << num
302  << " Nmod= " << numberOfModels << " " << this
303  << G4endl;
304  */
305  for(G4int i=0; i<numberOfModels; ++i) {
306  G4VMscModel* msc =
307  static_cast<G4VMscModel*>(GetModelByIndex(i, printing));
308  G4VMscModel* msc0=
309  static_cast<G4VMscModel*>(masterProcess->GetModelByIndex(i,printing));
310  msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false);
311  }
312  }
313 
314  }
315 
316  // explicitly defined printout by particle name
317  if(1 < verboseLevel ||
318  (0 < verboseLevel && (num == "e-" ||
319  num == "e+" || num == "mu+" ||
320  num == "mu-" || num == "proton"||
321  num == "pi+" || num == "pi-" ||
322  num == "kaon+" || num == "kaon-" ||
323  num == "alpha" || num == "anti_proton" ||
324  num == "GenericIon")))
325  {
326  G4cout << G4endl << GetProcessName()
327  << ": for " << num
328  << " SubType= " << GetProcessSubType()
329  << G4endl;
330  PrintInfo();
331  modelManager->DumpModelList(verboseLevel);
332  }
333 
334  if(1 < verboseLevel) {
335  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
336  << GetProcessName()
337  << " and particle " << num
338  << G4endl;
339  }
340 }
341 
342 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
343 
345 {
346  if (0 < verboseLevel) {
347  G4cout << G4endl << GetProcessName()
348  << ": for " << firstParticle->GetParticleName()
349  << " SubType= " << GetProcessSubType()
350  << G4endl;
351  PrintInfo();
352  modelManager->DumpModelList(verboseLevel);
353  }
354 }
355 
356 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
357 
359 {
360  G4VEnergyLossProcess* eloss = 0;
361  if(track->GetParticleDefinition() != currParticle) {
362  currParticle = track->GetParticleDefinition();
363  fIonisation = emManager->GetEnergyLossProcess(currParticle);
364  eloss = fIonisation;
365  }
366  /*
367  G4cout << "G4VMultipleScattering::StartTracking Nmod= " << numberOfModels
368  << " " << currParticle->GetParticleName()
369  << " E(MeV)= " << track->GetKineticEnergy()
370  << " Ion= " << eloss << " " << fIonisation << " IsMaster= "
371  << G4LossTableManager::Instance()->IsMaster()
372  << G4endl;
373  */
374  // one model
375  if(1 == numberOfModels) {
376  currentModel->StartTracking(track);
377  if(eloss) { currentModel->SetIonisation(fIonisation, currParticle); }
378 
379  // many models
380  } else {
381  for(G4int i=0; i<numberOfModels; ++i) {
382  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i,true));
383  msc->StartTracking(track);
384  if(eloss) { msc->SetIonisation(fIonisation, currParticle); }
385  }
386  }
387 }
388 
389 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
390 
392  const G4Track& track,
393  G4double,
394  G4double currentMinimalStep,
395  G4double&,
396  G4GPILSelection* selection)
397 {
398  // get Step limit proposed by the process
399  *selection = NotCandidateForSelection;
400  physStepLimit = gPathLength = tPathLength = currentMinimalStep;
401 
402  G4double ekin = track.GetKineticEnergy();
403  /*
404  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
405  << " " << currParticle->GetParticleName()
406  << " currMod " << currentModel
407  << G4endl;
408  */
409  // isIon flag is used only to select a model
410  if(isIon) {
411  ekin *= proton_mass_c2/track.GetParticleDefinition()->GetPDGMass();
412  }
413 
414  // select new model
415  if(1 < numberOfModels) {
416  currentModel = static_cast<G4VMscModel*>(
417  SelectModel(ekin,track.GetMaterialCutsCouple()->GetIndex()));
418  }
419  // step limit
420  if(currentModel->IsActive(ekin) && gPathLength >= geomMin
421  && ekin >= lowestKinEnergy) {
422  isActive = true;
423  tPathLength =
424  currentModel->ComputeTruePathLengthLimit(track, gPathLength);
425  if (tPathLength < physStepLimit) {
426  *selection = CandidateForSelection;
427  }
428  } else { isActive = false; }
429 
430 
431  //if(currParticle->GetPDGMass() > GeV)
432  /*
433  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
434  << " " << currParticle->GetParticleName()
435  << " gPathLength= " << gPathLength
436  << " tPathLength= " << tPathLength
437  << " currentMinimalStep= " << currentMinimalStep
438  << " isActive " << isActive << G4endl;
439  */
440  return gPathLength;
441 }
442 
443 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
444 
445 G4double
448 {
449  *condition = Forced;
450  //*condition = NotForced;
451  return DBL_MAX;
452 }
453 
454 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
455 
458 {
461  fNewPosition = step.GetPostStepPoint()->GetPosition();
462  fParticleChange.ProposePosition(fNewPosition);
463  fPositionChanged = false;
464 
465  G4double geomLength = step.GetStepLength();
466 
467  // very small step - no msc
468  if(!isActive) {
469  tPathLength = geomLength;
470 
471  // sample msc
472  } else {
473  G4double range =
474  currentModel->GetRange(currParticle,track.GetKineticEnergy(),
475  track.GetMaterialCutsCouple());
476 
477  tPathLength = currentModel->ComputeTrueStepLength(geomLength);
478 
479  // protection against wrong t->g->t conversion
480  /*
481  if(currParticle->GetPDGMass() > GeV)
482  G4cout << "G4VMsc::AlongStepDoIt: GeomLength= "
483  << geomLength
484  << " tPathLength= " << tPathLength
485  << " physStepLimit= " << physStepLimit
486  << " dr= " << range - trueLength
487  << " ekin= " << track.GetKineticEnergy() << G4endl;
488  */
489  if (tPathLength > physStepLimit) {
490  tPathLength = physStepLimit;
491  }
492 
493  // do not sample scattering at the last or at a small step
494  if(tPathLength + geomMin < range && tPathLength > geomMin) {
495 
496  G4double preSafety = step.GetPreStepPoint()->GetSafety();
497  G4double maxDisp = (tPathLength + geomLength)*0.5;
498  G4double postSafety= preSafety - maxDisp;
499  G4bool safetyRecomputed = false;
500  if(postSafety < maxDisp) {
501  safetyRecomputed = true;
502  postSafety = safetyHelper->ComputeSafety(fNewPosition,maxDisp);
503  }
504  G4ThreeVector displacement = currentModel->SampleScattering(
505  step.GetPostStepPoint()->GetMomentumDirection(), postSafety);
506 
507  G4double r2 = displacement.mag2();
508 
509  //G4cout << "R= " << sqrt(r2) << " postSafety= " << postSafety
510  // << G4endl;
511 
512  // make correction for displacement
513  if(r2 > 0.0) {
514 
515  fPositionChanged = true;
516  G4double fac = 1.0;
517 
518  // displaced point is definitely within the volume
519  if(r2 > postSafety*postSafety) {
520  G4double dispR = std::sqrt(r2);
521  if(!safetyRecomputed) {
522  postSafety = safetyHelper->ComputeSafety(fNewPosition, dispR);
523  }
524 
525  if(dispR > postSafety) {
526  fac = 0.99*postSafety/dispR;
527  }
528  }
529  // compute new endpoint of the Step
530  fNewPosition += fac*displacement;
531  //safetyHelper->ReLocateWithinVolume(fNewPosition);
532  }
533  }
534  }
536  //fParticleChange.ProposePosition(fNewPosition);
537  return &fParticleChange;
538 }
539 
540 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
541 
544 {
545  fParticleChange.Initialize(track);
546 
547  if(fPositionChanged) {
548  safetyHelper->ReLocateWithinVolume(fNewPosition);
549  fParticleChange.ProposePosition(fNewPosition);
550  }
551 
552  return &fParticleChange;
553 }
554 
555 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
556 
558  const G4Track& track,
559  G4double previousStepSize,
560  G4double currentMinimalStep,
561  G4double& currentSafety)
562 {
564  G4double x = AlongStepGetPhysicalInteractionLength(track,previousStepSize,
565  currentMinimalStep,
566  currentSafety,
567  &selection);
568  return x;
569 }
570 
571 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
572 
574  const G4Track& track,
575  G4double previousStepSize,
576  G4double currentMinimalStep,
577  G4double& currentSafety)
578 {
579  return GetContinuousStepLimit(track,previousStepSize,currentMinimalStep,
580  currentSafety);
581 }
582 
583 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
584 
587 {
588  *condition = Forced;
589  return DBL_MAX;
590 }
591 
592 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
593 
594 G4bool
596  const G4String& directory,
597  G4bool ascii)
598 {
599  G4bool yes = true;
600  if(part != firstParticle) { return yes; }
601  G4int nmod = modelManager->NumberOfModels();
602  static const G4String ss[4] = {"1","2","3","4"};
603  for(G4int i=0; i<nmod; ++i) {
604  G4VEmModel* msc = modelManager->GetModel(i);
605  yes = true;
606  G4PhysicsTable* table = msc->GetCrossSectionTable();
607  if (table) {
608  G4int j = std::min(i,3);
609  G4String name =
610  GetPhysicsTableFileName(part,directory,"LambdaMod"+ss[j],ascii);
611  yes = table->StorePhysicsTable(name,ascii);
612 
613  if ( yes ) {
614  if ( verboseLevel>0 ) {
615  G4cout << "Physics table are stored for "
616  << part->GetParticleName()
617  << " and process " << GetProcessName()
618  << " with a name <" << name << "> " << G4endl;
619  }
620  } else {
621  G4cout << "Fail to store Physics Table for "
622  << part->GetParticleName()
623  << " and process " << GetProcessName()
624  << " in the directory <" << directory
625  << "> " << G4endl;
626  }
627  }
628  }
629  return yes;
630 }
631 
632 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
633 
634 G4bool
636  const G4String&,
637  G4bool)
638 {
639  return true;
640 }
641 
642 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
643 
645 {
646  for(G4int i=0; i<numberOfModels; ++i) {
647  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i, true));
648  msc->SetIonisation(p, firstParticle);
649  }
650 }
651 
652 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
653 
G4double condition(const G4ErrorSymMatrix &m)
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii)
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
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:271
void InitialiseHelper()
static G4LossTableManager * Instance()
G4int verboseLevel
Definition: G4VProcess.hh:368
Int_t index
Definition: macro.C:9
G4double GetStepLength() const
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:592
void SetLateralDisplasmentFlag(G4bool val)
G4VEmModel * GetModel(G4int, G4bool ver=false)
void DeRegister(G4VEnergyLossProcess *p)
const char * p
Definition: xmltok.h:285
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4double GetMeanFreePath(const G4Track &track, G4double, G4ForceCondition *condition)
void SetLateralDisplasmentFlag(G4bool val)
Definition: G4VMscModel.hh:196
G4double PostStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4ForceCondition *condition)
void SetEmModel(G4VMscModel *, G4int index=1)
const XML_Char * name
Definition: expat.h:151
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:784
tuple x
Definition: test.py:50
void SetStepLimitType(G4MscStepLimitType)
Definition: G4VMscModel.hh:224
G4ProcessManager * GetProcessManager() const
int G4int
Definition: G4Types.hh:78
const G4String & GetPhysicsTableFileName(const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
Definition: G4VProcess.cc:186
const G4String & GetParticleName() const
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:683
G4VMultipleScattering(const G4String &name="msc", G4ProcessType type=fElectromagnetic)
G4StepPoint * GetPreStepPoint() const
G4MscStepLimitType StepLimitType() const
G4ParticleChangeForMSC fParticleChange
Char_t n[5]
virtual G4double ComputeTrueStepLength(G4double geomPathLength)
Definition: G4VMscModel.cc:159
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
G4GLOB_DLL std::ostream G4cout
virtual G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &stepLimit)
Definition: G4VMscModel.cc:145
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
const G4ParticleDefinition const G4Material *G4double range
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=0)
Definition: G4VEmModel.cc:387
const G4ThreeVector & GetPosition() const
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:288
bool G4bool
Definition: G4Types.hh:79
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
void SetRangeFactor(G4double)
Definition: G4VMscModel.hh:210
void ProposeTrueStepLength(G4double truePathLength)
TString part[npart]
Definition: Style.C:32
const G4ParticleDefinition * GetParticleDefinition() const
G4bool IsMaster() const
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
float proton_mass_c2
Definition: hepunit.py:275
const G4String & GetParticleType() const
void SetCrossSectionTable(G4PhysicsTable *, G4bool isLocal)
Definition: G4VEmModel.cc:395
void Register(G4VEnergyLossProcess *p)
Definition: G4Step.hh:76
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
tuple v
Definition: test.py:18
string pname
Definition: eplot.py:33
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:711
void DumpModelList(G4int verb)
G4bool LateralDisplasmentFlag() const
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:669
G4int size() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
#define fm
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:335
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)
G4StepPoint * GetPostStepPoint() const
void AddEmModel(G4int order, G4VEmModel *, const G4Region *region=0)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double mag2() const
G4double GetSafety() const
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
void SetGeomFactor(G4double)
Definition: G4VMscModel.hh:217
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)
virtual void PrintInfo()=0
double G4double
Definition: G4Types.hh:76
void SetSkin(G4double)
Definition: G4VMscModel.hh:203
void BuildPhysicsTable(const G4ParticleDefinition &)
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4ForceCondition
void SetRangeFactor(G4double val)
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection)
#define DBL_MAX
Definition: templates.hh:83
G4double MaxKinEnergy() const
G4VMscModel * EmModel(G4int index=1) const
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:718
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:138
void SetIonisation(G4VEnergyLossProcess *)