Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4VMultipleScattering Class Referenceabstract

#include <G4VMultipleScattering.hh>

+ Inheritance diagram for G4VMultipleScattering:

Public Member Functions

 G4VMultipleScattering (const G4String &name="msc", G4ProcessType type=fElectromagnetic)
 
 ~G4VMultipleScattering () override
 
void ProcessDescription (std::ostream &outFile) const override
 
virtual void InitialiseProcess (const G4ParticleDefinition *)=0
 
void StreamInfo (std::ostream &outFile, const G4ParticleDefinition &, G4bool rst=false) const
 
void PreparePhysicsTable (const G4ParticleDefinition &) override
 
void BuildPhysicsTable (const G4ParticleDefinition &) override
 
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
 
G4double ContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
 
 G4VMultipleScattering (G4VMultipleScattering &)=delete
 
G4VMultipleScatteringoperator= (const G4VMultipleScattering &right)=delete
 
G4VEmModelSelectModel (G4double kinEnergy, size_t idx)
 
void AddEmModel (G4int order, G4VMscModel *, const G4Region *region=nullptr)
 
void SetEmModel (G4VMscModel *, G4int idx=0)
 
G4VMscModelEmModel (size_t index=0) const
 
G4int NumberOfModels () const
 
G4VMscModelGetModelByIndex (G4int idx, G4bool ver=false) const
 
G4bool LateralDisplasmentFlag () const
 
G4double Skin () const
 
G4double RangeFactor () const
 
G4double GeomFactor () const
 
G4double PolarAngleLimit () const
 
G4bool UseBaseMaterial () 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 ()
 
G4VContinuousDiscreteProcessoperator= (const G4VContinuousDiscreteProcess &)=delete
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
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 ()
 
G4VProcessoperator= (const G4VProcess &)=delete
 
G4bool operator== (const G4VProcess &right) const
 
G4bool operator!= (const G4VProcess &right) const
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
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)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
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 const G4VProcessGetCreatorProcess () const
 
virtual void StartTracking (G4Track *)
 
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
 
virtual void ProcessDescription (std::ostream &outfile) 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 StreamProcessInfo (std::ostream &) const
 
G4double GetMeanFreePath (const G4Track &track, G4double, G4ForceCondition *condition) override
 
G4double GetContinuousStepLimit (const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety) override
 
- Protected Member Functions inherited from G4VContinuousDiscreteProcess
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
virtual G4double GetContinuousStepLimit (const G4Track &aTrack, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)=0
 
void SetGPILSelection (G4GPILSelection selection)
 
G4GPILSelection GetGPILSelection () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Protected Attributes

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

Additional Inherited Members

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

Detailed Description

Definition at line 86 of file G4VMultipleScattering.hh.

Constructor & Destructor Documentation

◆ G4VMultipleScattering() [1/2]

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

Definition at line 86 of file G4VMultipleScattering.cc.

88 fNewPosition(0.,0.,0.),
89 fNewDirection(0.,0.,1.)
90{
91 theParameters = G4EmParameters::Instance();
94
95 lowestKinEnergy = 10*CLHEP::eV;
96
97 geomMin = 0.05*CLHEP::nm;
98 minDisplacement2 = geomMin*geomMin;
99
101
102 modelManager = new G4EmModelManager();
103 emManager = G4LossTableManager::Instance();
104 mscModels.reserve(2);
105 emManager->Register(this);
106}
@ fMultipleScattering
@ fElectromagnetic
static G4EmParameters * Instance()
static G4LossTableManager * Instance()
void Register(G4VEnergyLossProcess *p)
G4ParticleChangeForMSC fParticleChange
void SetVerboseLevel(G4int value)
Definition: G4VProcess.hh:416
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:410
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:325

◆ ~G4VMultipleScattering()

G4VMultipleScattering::~G4VMultipleScattering ( )
override

Definition at line 110 of file G4VMultipleScattering.cc.

111{
112 delete modelManager;
113 emManager->DeRegister(this);
114}
void DeRegister(G4VEnergyLossProcess *p)

◆ G4VMultipleScattering() [2/2]

G4VMultipleScattering::G4VMultipleScattering ( G4VMultipleScattering )
delete

Member Function Documentation

◆ AddEmModel()

void G4VMultipleScattering::AddEmModel ( G4int  order,
G4VMscModel ptr,
const G4Region region = nullptr 
)

Definition at line 118 of file G4VMultipleScattering.cc.

120{
121 if(nullptr == ptr) { return; }
122 G4VEmFluctuationModel* fm = nullptr;
123 modelManager->AddEmModel(order, ptr, fm, region);
125}
void AddEmModel(G4int, G4VEmModel *, G4VEmFluctuationModel *fm, const G4Region *r)
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:390

Referenced by G4AdjointhMultipleScattering::InitialiseProcess(), G4eAdjointMultipleScattering::InitialiseProcess(), G4MuMultipleScattering::InitialiseProcess(), G4eMultipleScattering::InitialiseProcess(), G4hMultipleScattering::InitialiseProcess(), and G4EmConfigurator::PrepareModels().

◆ AlongStepDoIt()

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 290 of file G4VMultipleScattering.cc.

291{
292 fParticleChange.InitialiseMSC(track, step);
293 fNewPosition = fParticleChange.GetProposedPosition();
294 fPositionChanged = false;
295
296 G4double geomLength = step.GetStepLength();
297
298 // very small step - no msc
299 if(!isActive) {
300 tPathLength = geomLength;
301
302 // sample msc
303 } else {
304 G4double range =
305 currentModel->GetRange(currParticle,track.GetKineticEnergy(),
306 track.GetMaterialCutsCouple());
307
308 tPathLength = currentModel->ComputeTrueStepLength(geomLength);
309
310 /*
311 if(currParticle->GetPDGMass() > 0.9*GeV)
312 G4cout << "G4VMsc::AlongStepDoIt: GeomLength= "
313 << geomLength
314 << " tPathLength= " << tPathLength
315 << " physStepLimit= " << physStepLimit
316 << " dr= " << range - tPathLength
317 << " ekin= " << track.GetKineticEnergy() << G4endl;
318 */
319 // protection against wrong t->g->t conversion
320 tPathLength = std::min(tPathLength, physStepLimit);
321
322 // do not sample scattering at the last or at a small step
323 if(tPathLength < range && tPathLength > geomMin) {
324
325 static const G4double minSafety = 1.20*CLHEP::nm;
326 static const G4double sFact = 0.99;
327
328 G4ThreeVector displacement = currentModel->SampleScattering(
329 step.GetPostStepPoint()->GetMomentumDirection(),minSafety);
330
331 G4double r2 = displacement.mag2();
332 //G4cout << " R= " << sqrt(r2) << " Rmin= " << sqrt(minDisplacement2)
333 // << " flag= " << fDispBeyondSafety << G4endl;
334 if(r2 > minDisplacement2) {
335
336 fPositionChanged = true;
337 G4double dispR = std::sqrt(r2);
338 G4double postSafety =
339 sFact*safetyHelper->ComputeSafety(fNewPosition, dispR);
340 //G4cout<<" R= "<< dispR<<" postSafety= "<<postSafety<<G4endl;
341
342 // far away from geometry boundary
343 if(postSafety > 0.0 && dispR <= postSafety) {
344 fNewPosition += displacement;
345
346 //near the boundary
347 } else {
348 // displaced point is definitely within the volume
349 //G4cout<<" R= "<<dispR<<" postSafety= "<<postSafety<<G4endl;
350 if(dispR < postSafety) {
351 fNewPosition += displacement;
352
353 // reduced displacement
354 } else if(postSafety > geomMin) {
355 fNewPosition += displacement*(postSafety/dispR);
356
357 // very small postSafety
358 } else {
359 fPositionChanged = false;
360 }
361 }
362 if(fPositionChanged) {
363 safetyHelper->ReLocateWithinVolume(fNewPosition);
364 fParticleChange.ProposePosition(fNewPosition);
365 }
366 }
367 }
368 }
370 return &fParticleChange;
371}
double G4double
Definition: G4Types.hh:83
double mag2() const
const G4ThreeVector & GetProposedPosition() const
void InitialiseMSC(const G4Track &, const G4Step &step)
void ProposePosition(const G4ThreeVector &finalPosition)
G4double ComputeSafety(const G4ThreeVector &pGlobalPoint, G4double maxRadius=DBL_MAX)
void ReLocateWithinVolume(const G4ThreeVector &pGlobalPoint)
const G4ThreeVector & GetMomentumDirection() const
G4double GetStepLength() const
G4StepPoint * GetPostStepPoint() const
G4double GetKineticEnergy() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
virtual G4double ComputeTrueStepLength(G4double geomPathLength)=0
G4double GetRange(const G4ParticleDefinition *part, G4double kineticEnergy, const G4MaterialCutsCouple *couple)
Definition: G4VMscModel.cc:188
virtual G4ThreeVector & SampleScattering(const G4ThreeVector &, G4double safety)=0
void ProposeTrueStepLength(G4double truePathLength)

◆ AlongStepGetPhysicalInteractionLength()

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 219 of file G4VMultipleScattering.cc.

225{
226 // get Step limit proposed by the process
227 *selection = NotCandidateForSelection;
228 physStepLimit = gPathLength = tPathLength = currentMinimalStep;
229
230 G4double ekin = track.GetKineticEnergy();
231 /*
232 G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
233 << " " << currParticle->GetParticleName()
234 << " currMod " << currentModel
235 << G4endl;
236 */
237 // isIon flag is used only to select a model
238 if(isIon) {
239 ekin *= proton_mass_c2/track.GetParticleDefinition()->GetPDGMass();
240 }
241 const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple();
242
243 // select new model, static cast is possible in this class
244 if(1 < numberOfModels) {
245 currentModel =
246 static_cast<G4VMscModel*>(SelectModel(ekin,couple->GetIndex()));
247 }
248 currentModel->SetCurrentCouple(couple);
249 // msc is active is model is active, energy above the limit,
250 // and step size is above the limit;
251 // if it is active msc may limit the step
252 if(currentModel->IsActive(ekin) && tPathLength > geomMin
253 && ekin >= lowestKinEnergy) {
254 isActive = true;
255 tPathLength =
256 currentModel->ComputeTruePathLengthLimit(track, gPathLength);
257 if (tPathLength < physStepLimit) {
258 *selection = CandidateForSelection;
259 }
260 } else {
261 isActive = false;
262 gPathLength = DBL_MAX;
263 }
264
265 //if(currParticle->GetPDGMass() > GeV)
266 /*
267 G4cout << "MSC::AlongStepGPIL: Ekin= " << ekin
268 << " " << currParticle->GetParticleName()
269 << " gPathLength= " << gPathLength
270 << " tPathLength= " << tPathLength
271 << " currentMinimalStep= " << currentMinimalStep
272 << " isActive " << isActive << G4endl;
273 */
274 return gPathLength;
275}
@ CandidateForSelection
@ NotCandidateForSelection
const G4ParticleDefinition * GetParticleDefinition() const
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:468
G4bool IsActive(G4double kinEnergy) const
Definition: G4VEmModel.hh:774
virtual G4double ComputeTruePathLengthLimit(const G4Track &track, G4double &stepLimit)=0
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
#define DBL_MAX
Definition: templates.hh:62

Referenced by GetContinuousStepLimit().

◆ BuildPhysicsTable()

void G4VMultipleScattering::BuildPhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 170 of file G4VMultipleScattering.cc.

171{
172 G4bool master = emManager->IsMaster();
173
174 if(firstParticle == &part) {
175 emManager->BuildPhysicsTable(firstParticle);
176 }
177 const G4VMultipleScattering* ptr = this;
178 if(!master) {
179 ptr = static_cast<const G4VMultipleScattering*>(GetMasterProcess());
180 }
181
182 G4EmTableUtil::BuildMscProcess(this, ptr, part, firstParticle,
183 numberOfModels, master);
184}
bool G4bool
Definition: G4Types.hh:86
static void BuildMscProcess(G4VMultipleScattering *proc, const G4VMultipleScattering *masterProc, const G4ParticleDefinition &part, const G4ParticleDefinition *firstPart, G4int nModels, G4bool master)
void BuildPhysicsTable(const G4ParticleDefinition *aParticle)
const G4VProcess * GetMasterProcess() const
Definition: G4VProcess.hh:522

◆ ContinuousStepLimit()

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

Definition at line 391 of file G4VMultipleScattering.cc.

396{
397 return GetContinuousStepLimit(track,previousStepSize,currentMinimalStep,
398 currentSafety);
399}
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety) override

◆ EmModel()

G4VMscModel * G4VMultipleScattering::EmModel ( size_t  index = 0) const
inline

◆ FirstParticle()

const G4ParticleDefinition * G4VMultipleScattering::FirstParticle ( ) const
inline

Definition at line 361 of file G4VMultipleScattering.hh.

362{
363 return firstParticle;
364}

◆ GeomFactor()

G4double G4VMultipleScattering::GeomFactor ( ) const
inline

◆ GetContinuousStepLimit()

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

Implements G4VContinuousDiscreteProcess.

Definition at line 375 of file G4VMultipleScattering.cc.

380{
382 G4double x = AlongStepGetPhysicalInteractionLength(track,previousStepSize,
383 currentMinimalStep,
384 currentSafety,
385 &selection);
386 return x;
387}
G4GPILSelection
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection) override

Referenced by ContinuousStepLimit().

◆ GetMeanFreePath()

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

Implements G4VContinuousDiscreteProcess.

Definition at line 403 of file G4VMultipleScattering.cc.

405{
406 *condition = Forced;
407 return DBL_MAX;
408}
G4double condition(const G4ErrorSymMatrix &m)
@ Forced

◆ GetModelByIndex()

G4VMscModel * G4VMultipleScattering::GetModelByIndex ( G4int  idx,
G4bool  ver = false 
) const
inline

Definition at line 383 of file G4VMultipleScattering.hh.

384{
385 // static cast is possible inside this class
386 return static_cast<G4VMscModel*>(modelManager->GetModel(idx, ver));
387}
G4VEmModel * GetModel(G4int idx, G4bool ver=false) const

Referenced by G4EmTableUtil::BuildMscProcess(), G4EmTableUtil::PrepareMscProcess(), PreparePhysicsTable(), StartTracking(), and G4EmTableUtil::StoreMscTable().

◆ InitialiseProcess()

virtual void G4VMultipleScattering::InitialiseProcess ( const G4ParticleDefinition )
pure virtual

◆ LateralDisplasmentFlag()

◆ LowestKinEnergy()

G4double G4VMultipleScattering::LowestKinEnergy ( ) const
inline

Definition at line 347 of file G4VMultipleScattering.hh.

348{
349 return lowestKinEnergy;
350}

◆ NumberOfModels()

G4int G4VMultipleScattering::NumberOfModels ( ) const
inline

Definition at line 375 of file G4VMultipleScattering.hh.

376{
377 return numberOfModels;
378}

◆ operator=()

G4VMultipleScattering & G4VMultipleScattering::operator= ( const G4VMultipleScattering right)
delete

◆ PolarAngleLimit()

G4double G4VMultipleScattering::PolarAngleLimit ( ) const
inline

Definition at line 326 of file G4VMultipleScattering.hh.

327{
328 return theParameters->MscThetaLimit();
329}
G4double MscThetaLimit() const

Referenced by G4MuMultipleScattering::StreamProcessInfo().

◆ PostStepGetPhysicalInteractionLength()

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 280 of file G4VMultipleScattering.cc.

282{
284 return DBL_MAX;
285}
@ NotForced

◆ PreparePhysicsTable()

void G4VMultipleScattering::PreparePhysicsTable ( const G4ParticleDefinition part)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 141 of file G4VMultipleScattering.cc.

142{
143 G4bool master = emManager->IsMaster();
144 if(nullptr == firstParticle) { firstParticle = &part; }
145
146 emManager->PreparePhysicsTable(&part, this, master);
147 currParticle = nullptr;
148
149 if(firstParticle == &part) {
150 baseMat = emManager->GetTableBuilder()->GetBaseMaterialFlag();
151
152 G4EmTableUtil::PrepareMscProcess(this, part, modelManager,
153 stepLimit, facrange,
154 latDisplacement, master,
155 isIon, baseMat);
156
157 numberOfModels = modelManager->NumberOfModels();
158 currentModel = GetModelByIndex(0);
159
160 if(nullptr == safetyHelper) {
162 ->GetSafetyHelper();
163 safetyHelper->InitialiseHelper();
164 }
165 }
166}
G4int NumberOfModels() const
static void PrepareMscProcess(G4VMultipleScattering *proc, const G4ParticleDefinition &part, G4EmModelManager *modelManager, G4MscStepLimitType &stepLimit, G4double &facrange, G4bool &latDisplacement, G4bool &master, G4bool &isIon, G4bool &baseMat)
G4LossTableBuilder * GetTableBuilder()
void PreparePhysicsTable(const G4ParticleDefinition *aParticle, G4VEnergyLossProcess *p, G4bool theMaster)
void InitialiseHelper()
static G4TransportationManager * GetTransportationManager()
G4SafetyHelper * GetSafetyHelper() const
G4VMscModel * GetModelByIndex(G4int idx, G4bool ver=false) const

◆ ProcessDescription()

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

Reimplemented from G4VProcess.

Definition at line 437 of file G4VMultipleScattering.cc.

438{
439 if(nullptr != firstParticle) {
440 StreamInfo(outFile, *firstParticle, true);
441 }
442}
void StreamInfo(std::ostream &outFile, const G4ParticleDefinition &, G4bool rst=false) const

Referenced by G4MuMultipleScattering::ProcessDescription(), G4eMultipleScattering::ProcessDescription(), and G4hMultipleScattering::ProcessDescription().

◆ RangeFactor()

◆ RetrievePhysicsTable()

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

Reimplemented from G4VProcess.

Definition at line 428 of file G4VMultipleScattering.cc.

431{
432 return true;
433}

◆ SelectModel()

G4VEmModel * G4VMultipleScattering::SelectModel ( G4double  kinEnergy,
size_t  idx 
)
inline

Definition at line 291 of file G4VMultipleScattering.hh.

292{
293 return modelManager->SelectModel(kinEnergy, coupleIndex);
294}
G4VEmModel * SelectModel(G4double energy, std::size_t index)

Referenced by AlongStepGetPhysicalInteractionLength().

◆ SetEmModel()

void G4VMultipleScattering::SetEmModel ( G4VMscModel ptr,
G4int  idx = 0 
)

◆ SetLowestKinEnergy()

void G4VMultipleScattering::SetLowestKinEnergy ( G4double  val)
inline

Definition at line 354 of file G4VMultipleScattering.hh.

355{
356 lowestKinEnergy = val;
357}

◆ SetStepLimitType()

void G4VMultipleScattering::SetStepLimitType ( G4MscStepLimitType  val)
inline

Definition at line 340 of file G4VMultipleScattering.hh.

341{
342 theParameters->SetMscStepLimitType(val);
343}
void SetMscStepLimitType(G4MscStepLimitType val)

Referenced by LBE::ConstructEM().

◆ Skin()

G4double G4VMultipleScattering::Skin ( ) const
inline

◆ StartTracking()

void G4VMultipleScattering::StartTracking ( G4Track track)
overridevirtual

Reimplemented from G4VProcess.

Definition at line 200 of file G4VMultipleScattering.cc.

201{
202 G4VEnergyLossProcess* eloss = nullptr;
203 if(track->GetParticleDefinition() != currParticle) {
204 currParticle = track->GetParticleDefinition();
205 fIonisation = emManager->GetEnergyLossProcess(currParticle);
206 eloss = fIonisation;
207 }
208 for(G4int i=0; i<numberOfModels; ++i) {
210 msc->StartTracking(track);
211 if(nullptr != eloss) {
212 msc->SetIonisation(eloss, currParticle);
213 }
214 }
215}
int G4int
Definition: G4Types.hh:85
G4VEnergyLossProcess * GetEnergyLossProcess(const G4ParticleDefinition *)
virtual void StartTracking(G4Track *)
Definition: G4VEmModel.cc:214
void SetIonisation(G4VEnergyLossProcess *, const G4ParticleDefinition *part)
Definition: G4VMscModel.hh:315

Referenced by G4eAdjointMultipleScattering::StartTracking().

◆ StepLimitType()

◆ StorePhysicsTable()

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

Reimplemented from G4VProcess.

Definition at line 413 of file G4VMultipleScattering.cc.

416{
417 G4bool yes = true;
418 if(part != firstParticle || !emManager->IsMaster()) { return yes; }
419
420 return G4EmTableUtil::StoreMscTable(this, part, directory,
421 numberOfModels, verboseLevel,
422 ascii);
423}
static G4bool StoreMscTable(G4VMultipleScattering *proc, const G4ParticleDefinition *part, const G4String &directory, const G4int nModels, const G4int verb, const G4bool ascii)
G4int verboseLevel
Definition: G4VProcess.hh:360

◆ StreamInfo()

void G4VMultipleScattering::StreamInfo ( std::ostream &  outFile,
const G4ParticleDefinition part,
G4bool  rst = false 
) const

Definition at line 188 of file G4VMultipleScattering.cc.

190{
191 G4String indent = (rst ? " " : "");
192 outFile << G4endl << indent << GetProcessName() << ": ";
193 if (!rst) outFile << " for " << part.GetParticleName();
194 outFile << " SubType= " << GetProcessSubType() << G4endl;
195 modelManager->DumpModelList(outFile, verboseLevel);
196}
#define G4endl
Definition: G4ios.hh:57
void DumpModelList(std::ostream &out, G4int verb)
const G4String & GetParticleName() const
G4int GetProcessSubType() const
Definition: G4VProcess.hh:404
const G4String & GetProcessName() const
Definition: G4VProcess.hh:386

Referenced by G4EmTableUtil::BuildMscProcess(), and ProcessDescription().

◆ StreamProcessInfo()

virtual void G4VMultipleScattering::StreamProcessInfo ( std::ostream &  ) const
inlineprotectedvirtual

◆ UseBaseMaterial()

G4bool G4VMultipleScattering::UseBaseMaterial ( ) const
inline

Definition at line 391 of file G4VMultipleScattering.hh.

392{
393 return baseMat;
394}

Referenced by G4EmTableUtil::BuildMscProcess().

Member Data Documentation

◆ fParticleChange

G4ParticleChangeForMSC G4VMultipleScattering::fParticleChange
protected

Definition at line 260 of file G4VMultipleScattering.hh.

Referenced by AlongStepDoIt(), and G4VMultipleScattering().


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