Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VMultipleScattering Class Referenceabstract

#include <G4VMultipleScattering.hh>

Inheritance diagram for G4VMultipleScattering:
Collaboration diagram for G4VMultipleScattering:

Public Member Functions

 G4VMultipleScattering (const G4String &name="msc", G4ProcessType type=fElectromagnetic)
 
virtual ~G4VMultipleScattering ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &p) override=0
 
virtual void PrintInfo ()=0
 
virtual void ProcessDescription (std::ostream &outFile) const
 
void PreparePhysicsTable (const G4ParticleDefinition &) override
 
void BuildPhysicsTable (const G4ParticleDefinition &) override
 
void PrintInfoDefinition ()
 
G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
 
G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
 
void StartTracking (G4Track *) override
 
G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection) override
 
G4double PostStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4ForceCondition *condition) override
 
G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &) override
 
G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &) override
 
G4double ContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
 
G4VEmModelSelectModel (G4double kinEnergy, size_t idx)
 
void AddEmModel (G4int order, G4VEmModel *, const G4Region *region=nullptr)
 
void SetEmModel (G4VMscModel *, G4int index=1)
 
G4VMscModelEmModel (G4int index=1) const
 
G4VEmModelGetModelByIndex (G4int idx=0, G4bool ver=false) const
 
void SetIonisation (G4VEnergyLossProcess *)
 
G4bool LateralDisplasmentFlag () const
 
void SetLateralDisplasmentFlag (G4bool val)
 
G4double Skin () const
 
void SetSkin (G4double val)
 
G4double RangeFactor () const
 
void SetRangeFactor (G4double val)
 
G4double GeomFactor () const
 
G4double PolarAngleLimit () const
 
G4MscStepLimitType StepLimitType () const
 
void SetStepLimitType (G4MscStepLimitType val)
 
G4double LowestKinEnergy () const
 
void SetLowestKinEnergy (G4double val)
 
const G4ParticleDefinitionFirstParticle () const
 
- Public Member Functions inherited from G4VContinuousDiscreteProcess
 G4VContinuousDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VContinuousDiscreteProcess (G4VContinuousDiscreteProcess &)
 
virtual ~G4VContinuousDiscreteProcess ()
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual void InitialiseProcess (const G4ParticleDefinition *)=0
 
G4double GetMeanFreePath (const G4Track &track, G4double, G4ForceCondition *condition) override
 
G4double GetContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety) override
 
size_t NumberOfModels () const
 
- Protected Member Functions inherited from G4VContinuousDiscreteProcess
void SetGPILSelection (G4GPILSelection selection)
 
G4GPILSelection GetGPILSelection () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Protected Attributes

G4ParticleChangeForMSC fParticleChange
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 

Detailed Description

Definition at line 92 of file G4VMultipleScattering.hh.

Constructor & Destructor Documentation

G4VMultipleScattering::G4VMultipleScattering ( const G4String name = "msc",
G4ProcessType  type = fElectromagnetic 
)

Definition at line 91 of file G4VMultipleScattering.cc.

93  numberOfModels(0),
94  firstParticle(nullptr),
95  currParticle(nullptr),
96  stepLimit(fUseSafety),
97  facrange(0.04),
98  latDisplacement(true),
99  isIon(false),
100  fDispBeyondSafety(false),
101  fNewPosition(0.,0.,0.),
102  fNewDirection(0.,0.,1.)
103 {
104  theParameters = G4EmParameters::Instance();
105  SetVerboseLevel(1);
107  if("ionmsc" == name) { firstParticle = G4GenericIon::GenericIon(); }
108 
109  lowestKinEnergy = 10*CLHEP::eV;
110 
111  physStepLimit = gPathLength = tPathLength = 0.0;
112  fIonisation = nullptr;
113 
114  geomMin = 0.05*CLHEP::nm;
115  minDisplacement2 = geomMin*geomMin;
116 
118  safetyHelper = nullptr;
119  fPositionChanged = false;
120  isActive = false;
121 
122  currentModel = nullptr;
123  modelManager = new G4EmModelManager();
124  emManager = G4LossTableManager::Instance();
125  emManager->Register(this);
126 }
static G4LossTableManager * Instance()
G4ParticleChangeForMSC fParticleChange
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
void Register(G4VEnergyLossProcess *p)
static constexpr double eV
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
static constexpr double nm
Definition: SystemOfUnits.h:92
static G4EmParameters * Instance()
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437

Here is the call graph for this function:

G4VMultipleScattering::~G4VMultipleScattering ( )
virtual

Definition at line 130 of file G4VMultipleScattering.cc.

131 {
132  /*
133  if(1 < verboseLevel) {
134  G4cout << "G4VMultipleScattering destruct " << GetProcessName()
135  << G4endl;
136  }
137  */
138  delete modelManager;
139  emManager->DeRegister(this);
140 }
void DeRegister(G4VEnergyLossProcess *p)

Here is the call graph for this function:

Member Function Documentation

void G4VMultipleScattering::AddEmModel ( G4int  order,
G4VEmModel p,
const G4Region region = nullptr 
)

Definition at line 144 of file G4VMultipleScattering.cc.

146 {
147  G4VEmFluctuationModel* fm = nullptr;
148  modelManager->AddEmModel(order, p, fm, region);
149  if(p) { p->SetParticleChange(pParticleChange); }
150 }
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:418
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *, const G4Region *)
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283

Here is the call graph for this function:

Here is the caller graph for this function:

G4VParticleChange * G4VMultipleScattering::AlongStepDoIt ( const G4Track track,
const G4Step step 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 481 of file G4VMultipleScattering.cc.

482 {
485  fNewPosition = step.GetPostStepPoint()->GetPosition();
486  fParticleChange.ProposePosition(fNewPosition);
487  fPositionChanged = false;
488 
489  G4double geomLength = step.GetStepLength();
490 
491  // very small step - no msc
492  if(!isActive) {
493  tPathLength = geomLength;
494 
495  // sample msc
496  } else {
497  G4double range =
498  currentModel->GetRange(currParticle,track.GetKineticEnergy(),
499  track.GetMaterialCutsCouple());
500 
501  tPathLength = currentModel->ComputeTrueStepLength(geomLength);
502 
503  /*
504  if(currParticle->GetPDGMass() > 0.9*GeV)
505  G4cout << "G4VMsc::AlongStepDoIt: GeomLength= "
506  << geomLength
507  << " tPathLength= " << tPathLength
508  << " physStepLimit= " << physStepLimit
509  << " dr= " << range - tPathLength
510  << " ekin= " << track.GetKineticEnergy() << G4endl;
511  */
512  // protection against wrong t->g->t conversion
513  tPathLength = std::min(tPathLength, physStepLimit);
514 
515  // do not sample scattering at the last or at a small step
516  if(tPathLength < range && tPathLength > geomMin) {
517 
518  static const G4double minSafety = 1.20*CLHEP::nm;
519  static const G4double sFact = 0.99;
520 
521  G4ThreeVector displacement = currentModel->SampleScattering(
522  step.GetPostStepPoint()->GetMomentumDirection(),minSafety);
523 
524  G4double r2 = displacement.mag2();
525  //G4cout << " R= " << sqrt(r2) << " Rmin= " << sqrt(minDisplacement2)
526  // << " flag= " << fDispBeyondSafety << G4endl;
527  if(r2 > minDisplacement2) {
528 
529  fPositionChanged = true;
530  G4double dispR = std::sqrt(r2);
531  G4double postSafety =
532  sFact*safetyHelper->ComputeSafety(fNewPosition, dispR);
533  //G4cout<<" R= "<< dispR<<" postSafety= "<<postSafety<<G4endl;
534 
535  // far away from geometry boundary
536  if(postSafety > 0.0 && dispR <= postSafety) {
537  fNewPosition += displacement;
538 
539  //near the boundary
540  } else {
541  // displaced point is definitely within the volume
542  //G4cout<<" R= "<<dispR<<" postSafety= "<<postSafety<<G4endl;
543  if(dispR < postSafety) {
544  fNewPosition += displacement;
545 
546  // optional extra mechanism is applied only if a particle
547  // is stopped by the boundary
548  } else if(fDispBeyondSafety && 0.0 == postSafety) {
549  fNewPosition += displacement;
550  G4double maxshift =
551  std::min(2.0*dispR, geomLength*(physStepLimit/tPathLength - 1.0));
552  G4double dist = 0.0;
553  G4double safety = postSafety + dispR;
554  fNewDirection = *(fParticleChange.GetMomentumDirection());
555  /*
556  G4cout << "##MSC before Recheck maxshift= " << maxshift
557  << " postsafety= " << postSafety
558  << " Ekin= " << track.GetKineticEnergy()
559  << " " << track.GetDefinition()->GetParticleName()
560  << G4endl;
561  */
562  // check if it is possible to shift to the boundary
563  // and the shift is not large
564  if(safetyHelper->RecheckDistanceToCurrentBoundary(fNewPosition,
565  fNewDirection, maxshift, &dist, &safety)
566  && std::abs(dist) < maxshift) {
567  /*
568  G4cout << "##MSC after Recheck dist= " << dist
569  << " postsafety= " << postSafety
570  << " t= " << tPathLength
571  << " g= " << geomLength
572  << " p= " << physStepLimit
573  << G4endl;
574  */
575  // shift is positive
576  if(dist >= 0.0) {
577  tPathLength *= (1.0 + dist/geomLength);
578  fNewPosition += dist*fNewDirection;
579 
580  // shift is negative cannot be larger than geomLength
581  } else {
582  maxshift = std::min(maxshift, geomLength);
583  if(0.0 < maxshift + dist) {
584  G4ThreeVector postpoint = step.GetPostStepPoint()->GetPosition();
585  G4ThreeVector point = fNewPosition + dist*fNewDirection;
586  G4double R2 = (postpoint - point).mag2();
587  G4double newdist = dist;
588  // check not more than 10 extra boundaries
589  for(G4int i=0; i<10; ++i) {
590  dist = 0.0;
591  if(safetyHelper->RecheckDistanceToCurrentBoundary(
592  point, fNewDirection, maxshift, &dist, &safety)
593  && std::abs(newdist + dist) < maxshift) {
594  point += dist*fNewDirection;
595  G4double R2new = (postpoint - point).mag2();
596  //G4cout << "Backward i= " << i << " dist= " << dist
597  // << " R2= " << R2new << G4endl;
598  if(dist >= 0.0 || R2new > R2) { break; }
599  R2 = R2new;
600  fNewPosition = point;
601  newdist += dist;
602  } else {
603  break;
604  }
605  }
606  tPathLength *= (1.0 + newdist/geomLength);
607  // shift on boundary is not possible for negative disp
608  } else {
609  fNewPosition += displacement*(postSafety/dispR - 1.0);
610  }
611  }
612  // shift on boundary is not possible for any disp
613  } else {
614  fNewPosition += displacement*(postSafety/dispR - 1.0);
615  }
616  // reduced displacement
617  } else if(postSafety > geomMin) {
618  fNewPosition += displacement*(postSafety/dispR);
619 
620  // very small postSafety
621  } else {
622  fPositionChanged = false;
623  }
624  }
625  if(fPositionChanged) {
626  safetyHelper->ReLocateWithinVolume(fNewPosition);
627  fParticleChange.ProposePosition(fNewPosition);
628  }
629  }
630  }
631  }
633  return &fParticleChange;
634 }
const G4ThreeVector * GetMomentumDirection() const
void ReLocateWithinVolume(const G4ThreeVector &pGlobalPoint)
G4double GetStepLength() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
int G4int
Definition: G4Types.hh:78
G4ParticleChangeForMSC fParticleChange
virtual G4double ComputeTrueStepLength(G4double geomPathLength)
Definition: G4VMscModel.cc:163
const G4ThreeVector & GetMomentumDirection() const
G4double GetKineticEnergy() const
const G4ParticleDefinition const G4Material *G4double range
const G4ThreeVector & GetPosition() const
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.hh:283
void ProposeTrueStepLength(G4double truePathLength)
G4double ComputeSafety(const G4ThreeVector &pGlobalPoint, G4double maxRadius=DBL_MAX)
void ProposeMomentumDirection(const G4ThreeVector &Pfinal)
void ProposePosition(const G4ThreeVector &finalPosition)
static constexpr double nm
Definition: SystemOfUnits.h:92
G4StepPoint * GetPostStepPoint() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double mag2() const
double G4double
Definition: G4Types.hh:76
G4bool RecheckDistanceToCurrentBoundary(const G4ThreeVector &pGlobalPoint, const G4ThreeVector &pDirection, const G4double pCurrentProposedStepLength, G4double *prDistance, G4double *prNewSafety=0) const
virtual G4ThreeVector & SampleScattering(const G4ThreeVector &, G4double safety)
Definition: G4VMscModel.cc:142

Here is the call graph for this function:

G4double G4VMultipleScattering::AlongStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimalStep,
G4double currentSafety,
G4GPILSelection selection 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 414 of file G4VMultipleScattering.cc.

420 {
421  // get Step limit proposed by the process
422  *selection = NotCandidateForSelection;
423  physStepLimit = gPathLength = tPathLength = currentMinimalStep;
424 
425  G4double ekin = track.GetKineticEnergy();
426  /*
427  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
428  << " " << currParticle->GetParticleName()
429  << " currMod " << currentModel
430  << G4endl;
431  */
432  // isIon flag is used only to select a model
433  if(isIon) {
434  ekin *= proton_mass_c2/track.GetParticleDefinition()->GetPDGMass();
435  }
436 
437  // select new model
438  if(1 < numberOfModels) {
439  currentModel = static_cast<G4VMscModel*>(
440  SelectModel(ekin,track.GetMaterialCutsCouple()->GetIndex()));
441  }
442  // msc is active is model is active, energy above the limit,
443  // and step size is above the limit;
444  // if it is active msc may limit the step
445  if(currentModel->IsActive(ekin) && tPathLength > geomMin
446  && ekin >= lowestKinEnergy) {
447  isActive = true;
448  tPathLength =
449  currentModel->ComputeTruePathLengthLimit(track, gPathLength);
450  if (tPathLength < physStepLimit) {
451  *selection = CandidateForSelection;
452  }
453  } else { isActive = false; }
454 
455  //if(currParticle->GetPDGMass() > GeV)
456  /*
457  G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
458  << " " << currParticle->GetParticleName()
459  << " gPathLength= " << gPathLength
460  << " tPathLength= " << tPathLength
461  << " currentMinimalStep= " << currentMinimalStep
462  << " isActive " << isActive << G4endl;
463  */
464  return gPathLength;
465 }
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4double GetKineticEnergy() const
virtual G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &stepLimit)
Definition: G4VMscModel.cc:149
const G4ParticleDefinition * GetParticleDefinition() const
float proton_mass_c2
Definition: hepunit.py:275
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:760
G4double GetPDGMass() const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VMultipleScattering::BuildPhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 291 of file G4VMultipleScattering.cc.

292 {
293  G4String num = part.GetParticleName();
294  if(1 < verboseLevel) {
295  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
296  << GetProcessName()
297  << " and particle " << num
298  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
299  << G4endl;
300  }
301  G4bool master = true;
302  const G4VMultipleScattering* masterProcess =
303  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
304  if(masterProcess && masterProcess != this) { master = false; }
305 
306  if(firstParticle == &part) {
307  /*
308  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
309  << GetProcessName()
310  << " and particle " << num
311  << " IsMaster= " << G4LossTableManager::Instance()->IsMaster()
312  << " " << this
313  << G4endl;
314  */
315  emManager->BuildPhysicsTable(firstParticle);
316 
317  if(!master) {
318  // initialisation of models
319  G4bool printing = true;
320  numberOfModels = modelManager->NumberOfModels();
321  /*
322  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() for "
323  << GetProcessName()
324  << " and particle " << num
325  << " Nmod= " << numberOfModels << " " << this
326  << G4endl;
327  */
328  for(G4int i=0; i<numberOfModels; ++i) {
329  G4VMscModel* msc =
330  static_cast<G4VMscModel*>(GetModelByIndex(i, printing));
331  G4VMscModel* msc0=
332  static_cast<G4VMscModel*>(masterProcess->GetModelByIndex(i,printing));
333  msc->SetCrossSectionTable(msc0->GetCrossSectionTable(), false);
334  msc->InitialiseLocal(firstParticle, msc0);
335  }
336  }
337  }
338 
339  // explicitly defined printout by particle name
340  if(1 < verboseLevel ||
341  (0 < verboseLevel && (num == "e-" ||
342  num == "e+" || num == "mu+" ||
343  num == "mu-" || num == "proton"||
344  num == "pi+" || num == "pi-" ||
345  num == "kaon+" || num == "kaon-" ||
346  num == "alpha" || num == "anti_proton" ||
347  num == "GenericIon")))
348  {
349  G4cout << G4endl << GetProcessName()
350  << ": for " << num
351  << " SubType= " << GetProcessSubType()
352  << G4endl;
353  PrintInfo();
354  modelManager->DumpModelList(verboseLevel);
355  }
356 
357  if(1 < verboseLevel) {
358  G4cout << "### G4VMultipleScattering::BuildPhysicsTable() done for "
359  << GetProcessName()
360  << " and particle " << num
361  << G4endl;
362  }
363 }
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
static G4LossTableManager * Instance()
G4int verboseLevel
Definition: G4VProcess.hh:368
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:218
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:833
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
G4bool IsMaster() const
void SetCrossSectionTable(G4PhysicsTable *, G4bool isLocal)
Definition: G4VEmModel.cc:426
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void DumpModelList(G4int verb)
G4int NumberOfModels() const
void BuildPhysicsTable(const G4ParticleDefinition *aParticle)
#define G4endl
Definition: G4ios.hh:61
virtual void PrintInfo()=0
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const

Here is the call graph for this function:

G4double G4VMultipleScattering::ContinuousStepLimit ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimalStep,
G4double currentSafety 
)

Definition at line 663 of file G4VMultipleScattering.cc.

668 {
669  return GetContinuousStepLimit(track,previousStepSize,currentMinimalStep,
670  currentSafety);
671 }
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety) override

Here is the call graph for this function:

G4VMscModel * G4VMultipleScattering::EmModel ( G4int  index = 1) const

Definition at line 163 of file G4VMultipleScattering.cc.

164 {
165  G4VMscModel* p = nullptr;
166  if(index >= 0 && index < G4int(mscModels.size())) { p = mscModels[index]; }
167  return p;
168 }
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78

Here is the caller graph for this function:

const G4ParticleDefinition * G4VMultipleScattering::FirstParticle ( ) const
inline

Definition at line 405 of file G4VMultipleScattering.hh.

406 {
407  return firstParticle;
408 }
G4double G4VMultipleScattering::GeomFactor ( ) const
inline

Definition at line 361 of file G4VMultipleScattering.hh.

362 {
363  return theParameters->MscGeomFactor();
364 }
G4double MscGeomFactor() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VMultipleScattering::GetContinuousStepLimit ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimalStep,
G4double currentSafety 
)
overrideprotectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 647 of file G4VMultipleScattering.cc.

652 {
654  G4double x = AlongStepGetPhysicalInteractionLength(track,previousStepSize,
655  currentMinimalStep,
656  currentSafety,
657  &selection);
658  return x;
659 }
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection) override
tuple x
Definition: test.py:50
double G4double
Definition: G4Types.hh:76
G4GPILSelection

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VMultipleScattering::GetMeanFreePath ( const G4Track track,
G4double  ,
G4ForceCondition condition 
)
overrideprotectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 675 of file G4VMultipleScattering.cc.

677 {
678  *condition = Forced;
679  return DBL_MAX;
680 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83
G4VEmModel * G4VMultipleScattering::GetModelByIndex ( G4int  idx = 0,
G4bool  ver = false 
) const

Definition at line 173 of file G4VMultipleScattering.cc.

174 {
175  return modelManager->GetModel(idx, ver);
176 }
G4VEmModel * GetModel(G4int idx, G4bool ver=false)

Here is the call graph for this function:

Here is the caller graph for this function:

virtual void G4VMultipleScattering::InitialiseProcess ( const G4ParticleDefinition )
protectedpure virtual
virtual G4bool G4VMultipleScattering::IsApplicable ( const G4ParticleDefinition p)
overridepure virtual
G4bool G4VMultipleScattering::LateralDisplasmentFlag ( ) const
inline

Definition at line 316 of file G4VMultipleScattering.hh.

317 {
318  return latDisplacement;
319 }

Here is the caller graph for this function:

G4double G4VMultipleScattering::LowestKinEnergy ( ) const
inline

Definition at line 391 of file G4VMultipleScattering.hh.

392 {
393  return lowestKinEnergy;
394 }
size_t G4VMultipleScattering::NumberOfModels ( ) const
inlineprotected

Definition at line 412 of file G4VMultipleScattering.hh.

413 {
414  return modelManager->NumberOfModels();
415 }
G4int NumberOfModels() const

Here is the call graph for this function:

G4double G4VMultipleScattering::PolarAngleLimit ( ) const
inline

Definition at line 368 of file G4VMultipleScattering.hh.

369 {
370  return theParameters->MscThetaLimit();
371 }
G4double MscThetaLimit() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4VParticleChange * G4VMultipleScattering::PostStepDoIt ( const G4Track track,
const G4Step  
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 639 of file G4VMultipleScattering.cc.

640 {
642  return &fParticleChange;
643 }
virtual void Initialize(const G4Track &)
G4ParticleChangeForMSC fParticleChange

Here is the call graph for this function:

G4double G4VMultipleScattering::PostStepGetPhysicalInteractionLength ( const G4Track ,
G4double  previousStepSize,
G4ForceCondition condition 
)
overridevirtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 470 of file G4VMultipleScattering.cc.

472 {
473  //*condition = Forced;
474  *condition = NotForced;
475  return DBL_MAX;
476 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83
void G4VMultipleScattering::PreparePhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 181 of file G4VMultipleScattering.cc.

182 {
183  G4bool master = true;
184  const G4VMultipleScattering* masterProc =
185  static_cast<const G4VMultipleScattering*>(GetMasterProcess());
186  if(masterProc && masterProc != this) { master = false; }
187 
188  if(!firstParticle) { firstParticle = &part; }
189  if(part.GetParticleType() == "nucleus") {
190  stepLimit = fMinimal;
191  latDisplacement = false;
192  facrange = 0.2;
193  G4String pname = part.GetParticleName();
194  if(pname != "deuteron" && pname != "triton" &&
195  pname != "alpha+" && pname != "helium" &&
196  pname != "alpha" && pname != "He3" &&
197  pname != "hydrogen") {
198 
199  const G4ParticleDefinition* theGenericIon =
201  if(&part == theGenericIon) { isIon = true; }
202 
203  if(theGenericIon && firstParticle != theGenericIon) {
204  G4ProcessManager* pm = theGenericIon->GetProcessManager();
206  size_t n = v->size();
207  for(size_t j=0; j<n; ++j) {
208  if((*v)[j] == this) {
209  firstParticle = theGenericIon;
210  isIon = true;
211  break;
212  }
213  }
214  }
215  }
216  }
217 
218  emManager->PreparePhysicsTable(&part, this, master);
219  currParticle = 0;
220 
221  if(1 < verboseLevel) {
222  G4cout << "### G4VMultipleScattering::PrepearPhysicsTable() for "
223  << GetProcessName()
224  << " and particle " << part.GetParticleName()
225  << " local particle " << firstParticle->GetParticleName()
226  << " isIon= " << isIon
227  << G4endl;
228  }
229 
230  if(firstParticle == &part) {
231 
232  // initialise process
233  InitialiseProcess(firstParticle);
234 
235  // heavy particles and not ions
236  if(!isIon) {
237  if(part.GetPDGMass() > MeV) {
238  stepLimit = theParameters->MscMuHadStepLimitType();
239  facrange = theParameters->MscMuHadRangeFactor();
240  latDisplacement = theParameters->MuHadLateralDisplacement();
241  } else {
242  stepLimit = theParameters->MscStepLimitType();
243  facrange = theParameters->MscRangeFactor();
244  latDisplacement = theParameters->LateralDisplacement();
245  }
246  if(latDisplacement) {
247  fDispBeyondSafety = theParameters->LatDisplacementBeyondSafety();
248  }
249  }
250  if(master) { SetVerboseLevel(theParameters->Verbose()); }
251  else { SetVerboseLevel(theParameters->WorkerVerbose()); }
252 
253  // initialisation of models
254  numberOfModels = modelManager->NumberOfModels();
255  /*
256  G4cout << "### G4VMultipleScattering::PreparePhysicsTable() for "
257  << GetProcessName()
258  << " and particle " << part.GetParticleName()
259  << " Nmod= " << numberOfModels << " " << this
260  << G4endl;
261  */
262  for(G4int i=0; i<numberOfModels; ++i) {
263  G4VMscModel* msc = static_cast<G4VMscModel*>(modelManager->GetModel(i));
264  msc->SetIonisation(0, firstParticle);
265  msc->SetMasterThread(master);
266  if(0 == i) { currentModel = msc; }
267  msc->SetStepLimitType(stepLimit);
268  msc->SetLateralDisplasmentFlag(latDisplacement);
269  msc->SetSkin(theParameters->MscSkin());
270  msc->SetRangeFactor(facrange);
271  msc->SetGeomFactor(theParameters->MscGeomFactor());
272  msc->SetPolarAngleLimit(theParameters->MscThetaLimit());
273  G4double emax =
274  std::min(msc->HighEnergyLimit(),theParameters->MaxKinEnergy());
275  msc->SetHighEnergyLimit(emax);
276  }
277 
278  modelManager->Initialise(firstParticle, G4Electron::Electron(),
279  10.0, verboseLevel);
280 
281  if(!safetyHelper) {
283  ->GetSafetyHelper();
284  safetyHelper->InitialiseHelper();
285  }
286  }
287 }
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
G4SafetyHelper * GetSafetyHelper() const
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
void InitialiseHelper()
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VEmModel * GetModel(G4int idx, G4bool ver=false)
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:636
G4MscStepLimitType MscMuHadStepLimitType() const
G4double MscGeomFactor() const
G4double MscMuHadRangeFactor() const
G4double MscThetaLimit() const
void SetLateralDisplasmentFlag(G4bool val)
Definition: G4VMscModel.hh:195
void SetStepLimitType(G4MscStepLimitType)
Definition: G4VMscModel.hh:223
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4bool LatDisplacementBeyondSafety() const
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4GLOB_DLL std::ostream G4cout
G4int Verbose() const
G4bool MuHadLateralDisplacement() const
bool G4bool
Definition: G4Types.hh:79
void SetRangeFactor(G4double)
Definition: G4VMscModel.hh:209
G4bool LateralDisplacement() const
const G4String & GetParticleType() const
const G4int n
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
tuple v
Definition: test.py:18
string pname
Definition: eplot.py:33
const G4DataVector * Initialise(const G4ParticleDefinition *part, const G4ParticleDefinition *secPart, G4double minSubRange, G4int verb)
static G4TransportationManager * GetTransportationManager()
static const G4double emax
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:718
G4int size() const
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:330
G4int NumberOfModels() const
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
G4ProcessManager * GetProcessManager() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double MscRangeFactor() const
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
void SetGeomFactor(G4double)
Definition: G4VMscModel.hh:216
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4MscStepLimitType MscStepLimitType() const
double G4double
Definition: G4Types.hh:76
void SetSkin(G4double)
Definition: G4VMscModel.hh:202
G4ProcessVector * GetAlongStepProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
G4double MscSkin() const
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:767
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:437

Here is the call graph for this function:

virtual void G4VMultipleScattering::PrintInfo ( )
pure virtual
void G4VMultipleScattering::PrintInfoDefinition ( )

Definition at line 367 of file G4VMultipleScattering.cc.

368 {
369  if (0 < verboseLevel) {
370  G4cout << G4endl << GetProcessName()
371  << ": for " << firstParticle->GetParticleName()
372  << " SubType= " << GetProcessSubType()
373  << G4endl;
374  PrintInfo();
375  modelManager->DumpModelList(verboseLevel);
376  }
377 }
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void DumpModelList(G4int verb)
#define G4endl
Definition: G4ios.hh:61
virtual void PrintInfo()=0
G4int GetProcessSubType() const
Definition: G4VProcess.hh:426

Here is the call graph for this function:

void G4VMultipleScattering::ProcessDescription ( std::ostream &  outFile) const
virtual

Definition at line 748 of file G4VMultipleScattering.cc.

749 {
750  outFile << "Multiple scattering process <" << GetProcessName()
751  << ">" << G4endl;
752 }
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4double G4VMultipleScattering::RangeFactor ( ) const
inline

Definition at line 344 of file G4VMultipleScattering.hh.

345 {
346  return facrange;
347 }

Here is the caller graph for this function:

G4bool G4VMultipleScattering::RetrievePhysicsTable ( const G4ParticleDefinition ,
const G4String directory,
G4bool  ascii 
)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 729 of file G4VMultipleScattering.cc.

732 {
733  return true;
734 }
G4VEmModel * G4VMultipleScattering::SelectModel ( G4double  kinEnergy,
size_t  idx 
)
inline

Definition at line 309 of file G4VMultipleScattering.hh.

310 {
311  return modelManager->SelectModel(kinEnergy, coupleIndex);
312 }
G4VEmModel * SelectModel(G4double &energy, size_t &index)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VMultipleScattering::SetEmModel ( G4VMscModel p,
G4int  index = 1 
)

Definition at line 154 of file G4VMultipleScattering.cc.

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 }
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
const G4int n

Here is the caller graph for this function:

void G4VMultipleScattering::SetIonisation ( G4VEnergyLossProcess p)

Definition at line 738 of file G4VMultipleScattering.cc.

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 }
int G4int
Definition: G4Types.hh:78
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:330
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const

Here is the call graph for this function:

void G4VMultipleScattering::SetLateralDisplasmentFlag ( G4bool  val)
inline

Definition at line 323 of file G4VMultipleScattering.hh.

324 {
325  latDisplacement = val;
326 }
void G4VMultipleScattering::SetLowestKinEnergy ( G4double  val)
inline

Definition at line 398 of file G4VMultipleScattering.hh.

399 {
400  lowestKinEnergy = val;
401 }
void G4VMultipleScattering::SetRangeFactor ( G4double  val)
inline

Definition at line 351 of file G4VMultipleScattering.hh.

352 {
353  if(val > 0.0 && val < 1.0) {
354  facrange = val;
355  theParameters->SetMscRangeFactor(val);
356  }
357 }
void SetMscRangeFactor(G4double val)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VMultipleScattering::SetSkin ( G4double  val)
inline

Definition at line 337 of file G4VMultipleScattering.hh.

338 {
339  theParameters->SetMscSkin(val);
340 }
void SetMscSkin(G4double val)

Here is the call graph for this function:

void G4VMultipleScattering::SetStepLimitType ( G4MscStepLimitType  val)
inline

Definition at line 382 of file G4VMultipleScattering.hh.

383 {
384  stepLimit = val;
385  if(val == fMinimal) { SetRangeFactor(0.2); }
386  theParameters->SetMscStepLimitType(val);
387 }
void SetMscStepLimitType(G4MscStepLimitType val)
void SetRangeFactor(G4double val)

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4VMultipleScattering::Skin ( ) const
inline

Definition at line 330 of file G4VMultipleScattering.hh.

331 {
332  return theParameters->MscSkin();
333 }
G4double MscSkin() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4VMultipleScattering::StartTracking ( G4Track track)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 381 of file G4VMultipleScattering.cc.

382 {
383  G4VEnergyLossProcess* eloss = nullptr;
384  if(track->GetParticleDefinition() != currParticle) {
385  currParticle = track->GetParticleDefinition();
386  fIonisation = emManager->GetEnergyLossProcess(currParticle);
387  eloss = fIonisation;
388  }
389  /*
390  G4cout << "G4VMultipleScattering::StartTracking Nmod= " << numberOfModels
391  << " " << currParticle->GetParticleName()
392  << " E(MeV)= " << track->GetKineticEnergy()
393  << " Ion= " << eloss << " " << fIonisation << " IsMaster= "
394  << G4LossTableManager::Instance()->IsMaster()
395  << G4endl;
396  */
397  // one model
398  if(1 == numberOfModels) {
399  currentModel->StartTracking(track);
400  if(eloss) { currentModel->SetIonisation(fIonisation, currParticle); }
401 
402  // many models
403  } else {
404  for(G4int i=0; i<numberOfModels; ++i) {
405  G4VMscModel* msc = static_cast<G4VMscModel*>(GetModelByIndex(i,true));
406  msc->StartTracking(track);
407  if(eloss) { msc->SetIonisation(fIonisation, currParticle); }
408  }
409  }
410 }
virtual void StartTracking(G4Track *)
Definition: G4VEmModel.cc:292
int G4int
Definition: G4Types.hh:78
const G4ParticleDefinition * GetParticleDefinition() const
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:330
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const

Here is the call graph for this function:

Here is the caller graph for this function:

G4MscStepLimitType G4VMultipleScattering::StepLimitType ( ) const
inline

Definition at line 375 of file G4VMultipleScattering.hh.

376 {
377  return stepLimit;
378 }

Here is the caller graph for this function:

G4bool G4VMultipleScattering::StorePhysicsTable ( const G4ParticleDefinition part,
const G4String directory,
G4bool  ascii = false 
)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 685 of file G4VMultipleScattering.cc.

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 
695  G4int nmod = modelManager->NumberOfModels();
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 }
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:538
const XML_Char * name
Definition: expat.h:151
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VEmModel * GetModel(G4int idx, G4bool ver=false)
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:833
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
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4int NumberOfModels() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
#define G4endl
Definition: G4ios.hh:61
G4bool StorePhysicsTable(const G4String &filename, G4bool ascii=false)

Here is the call graph for this function:

Member Data Documentation

G4ParticleChangeForMSC G4VMultipleScattering::fParticleChange
protected

Definition at line 284 of file G4VMultipleScattering.hh.


The documentation for this class was generated from the following files: