Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4HadronElasticProcess Class Reference

#include <G4HadronElasticProcess.hh>

+ Inheritance diagram for G4HadronElasticProcess:

Public Member Functions

 G4HadronElasticProcess (const G4String &procName="hadElastic")
 
 ~G4HadronElasticProcess () override
 
G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep) override
 
virtual void SetLowestEnergy (G4double)
 
virtual void SetLowestEnergyNeutron (G4double)
 
void ProcessDescription (std::ostream &outFile) const override
 
void SetDiffraction (G4HadronicInteraction *, G4VCrossSectionRatio *)
 
- Public Member Functions inherited from G4HadronicProcess
 G4HadronicProcess (const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
 
 G4HadronicProcess (const G4String &processName, G4HadronicProcessType subType)
 
 ~G4HadronicProcess () override
 
void RegisterMe (G4HadronicInteraction *a)
 
G4double GetElementCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
G4double GetMicroscopicCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep) override
 
void PreparePhysicsTable (const G4ParticleDefinition &) override
 
void BuildPhysicsTable (const G4ParticleDefinition &) override
 
void DumpPhysicsTable (const G4ParticleDefinition &p)
 
void AddDataSet (G4VCrossSectionDataSet *aDataSet)
 
std::vector< G4HadronicInteraction * > & GetHadronicInteractionList ()
 
G4HadronicInteractionGetHadronicModel (const G4String &)
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *) override
 
const G4NucleusGetTargetNucleus () const
 
const G4IsotopeGetTargetIsotope ()
 
void ProcessDescription (std::ostream &outFile) const override
 
void BiasCrossSectionByFactor (G4double aScale)
 
void MultiplyCrossSectionBy (G4double factor)
 
G4double CrossSectionFactor () const
 
void SetIntegral (G4bool val)
 
void SetEpReportLevel (G4int level)
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
G4CrossSectionDataStoreGetCrossSectionDataStore ()
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &aName, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
G4VDiscreteProcessoperator= (const G4VDiscreteProcess &)=delete
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (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 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 &)
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Member Functions inherited from G4HadronicProcess
G4HadronicInteractionChooseHadronicInteraction (const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, const G4Material *aMaterial, const G4Element *anElement)
 
G4NucleusGetTargetNucleusPointer ()
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
G4HadFinalStateCheckResult (const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
 
void CheckEnergyMomentumConservation (const G4Track &, const G4Nucleus &)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4HadronicProcess
G4HadProjectile thePro
 
G4ParticleChangetheTotalResult
 
G4double fWeight
 
G4int epReportLevel
 
- 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
 

Detailed Description

Definition at line 49 of file G4HadronElasticProcess.hh.

Constructor & Destructor Documentation

◆ G4HadronElasticProcess()

G4HadronElasticProcess::G4HadronElasticProcess ( const G4String procName = "hadElastic")
explicit

Definition at line 46 of file G4HadronElasticProcess.cc.

48 fDiffraction(nullptr), fDiffractionRatio(nullptr)
49{}
@ fHadronElastic

◆ ~G4HadronElasticProcess()

G4HadronElasticProcess::~G4HadronElasticProcess ( )
override

Definition at line 51 of file G4HadronElasticProcess.cc.

52{}

Member Function Documentation

◆ PostStepDoIt()

G4VParticleChange * G4HadronElasticProcess::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
overridevirtual

Reimplemented from G4HadronicProcess.

Definition at line 62 of file G4HadronElasticProcess.cc.

64{
67 G4double weight = track.GetWeight();
69
70 // For elastic scattering, _any_ result is considered an interaction
72
73 G4double kineticEnergy = track.GetKineticEnergy();
74 G4TrackStatus status = track.GetTrackStatus();
75 if(kineticEnergy == 0.0 || track.GetTrackStatus() != fAlive) {
76 return theTotalResult;
77 }
78
79 const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
80 const G4ParticleDefinition* part = dynParticle->GetDefinition();
81 const G4Material* material = track.GetMaterial();
82 G4Nucleus* targNucleus = GetTargetNucleusPointer();
83
84 // Select element
85 const G4Element* elm =
86 GetCrossSectionDataStore()->SampleZandA(dynParticle, material, *targNucleus);
87
88 // Initialize the hadronic projectile from the track
89 G4HadProjectile theProj(track);
90 G4HadronicInteraction* hadi = nullptr;
91 G4HadFinalState* result = nullptr;
92
93 if(fDiffraction)
94 {
95 G4double ratio =
96 fDiffractionRatio->ComputeRatio(part, kineticEnergy,
97 targNucleus->GetZ_asInt(),
98 targNucleus->GetA_asInt());
99 // diffraction is chosen
100 if(ratio > 0.0 && G4UniformRand() < ratio)
101 {
102 try
103 {
104 result = fDiffraction->ApplyYourself(theProj, *targNucleus);
105 }
106 catch(G4HadronicException & aR)
107 {
109 aR.Report(ed);
110 ed << "Call for " << fDiffraction->GetModelName() << G4endl;
111 ed << "Target element "<< elm->GetName()<<" Z= "
112 << targNucleus->GetZ_asInt()
113 << " A= " << targNucleus->GetA_asInt() << G4endl;
114 DumpState(track,"ApplyYourself",ed);
115 ed << " ApplyYourself failed" << G4endl;
116 G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
117 FatalException, ed);
118 }
119 // Check the result for catastrophic energy non-conservation
120 result = CheckResult(theProj, *targNucleus, result);
121
122 result->SetTrafoToLab(theProj.GetTrafoToLab());
124
125 FillResult(result, track);
126
127 if (epReportLevel != 0) {
128 CheckEnergyMomentumConservation(track, *targNucleus);
129 }
130 return theTotalResult;
131 }
132 }
133
134 // ordinary elastic scattering
135 try
136 {
137 hadi = ChooseHadronicInteraction( theProj, *targNucleus, material, elm );
138 }
139 catch(G4HadronicException & aE)
140 {
142 aE.Report(ed);
143 ed << "Target element "<< elm->GetName()<<" Z= "
144 << targNucleus->GetZ_asInt() << " A= "
145 << targNucleus->GetA_asInt() << G4endl;
146 DumpState(track,"ChooseHadronicInteraction",ed);
147 ed << " No HadronicInteraction found out" << G4endl;
148 G4Exception("G4HadronElasticProcess::PostStepDoIt", "had005",
149 FatalException, ed);
150 }
151
152 size_t idx = track.GetMaterialCutsCouple()->GetIndex();
154 ->GetEnergyCutsVector(3)))[idx];
155 hadi->SetRecoilEnergyThreshold(tcut);
156
157 if(verboseLevel>1) {
158 G4cout << "G4HadronElasticProcess::PostStepDoIt for "
159 << part->GetParticleName()
160 << " in " << material->GetName()
161 << " Target Z= " << targNucleus->GetZ_asInt()
162 << " A= " << targNucleus->GetA_asInt()
163 << " Tcut(MeV)= " << tcut << G4endl;
164 }
165
166 try
167 {
168 result = hadi->ApplyYourself( theProj, *targNucleus);
169 }
170 catch(G4HadronicException & aR)
171 {
173 aR.Report(ed);
174 ed << "Call for " << hadi->GetModelName() << G4endl;
175 ed << "Target element "<< elm->GetName()<<" Z= "
176 << targNucleus->GetZ_asInt()
177 << " A= " << targNucleus->GetA_asInt() << G4endl;
178 DumpState(track,"ApplyYourself",ed);
179 ed << " ApplyYourself failed" << G4endl;
180 G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
181 FatalException, ed);
182 }
183
184 // Check the result for catastrophic energy non-conservation
185 // cannot be applied because is not guranteed that recoil
186 // nucleus is created
187 // result = CheckResult(theProj, targNucleus, result);
188
189 // directions
190 G4ThreeVector indir = track.GetMomentumDirection();
191 G4ThreeVector outdir = result->GetMomentumChange();
192
193 if(verboseLevel>1) {
194 G4cout << "Efin= " << result->GetEnergyChange()
195 << " de= " << result->GetLocalEnergyDeposit()
196 << " nsec= " << result->GetNumberOfSecondaries()
197 << " dir= " << outdir
198 << G4endl;
199 }
200
201 // energies
202 G4double edep = std::max(result->GetLocalEnergyDeposit(), 0.0);
203 G4double efinal = std::max(result->GetEnergyChange(), 0.0);
204
205 // primary change
207
208 if(efinal > 0.0) {
209 outdir.rotateUz(indir);
211 } else {
212 if(part->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
213 { status = fStopButAlive; }
214 else { status = fStopAndKill; }
216 }
217
218 //G4cout << "Efinal= " << efinal << " TrackStatus= " << status << G4endl;
220
221 // recoil
222 if(result->GetNumberOfSecondaries() > 0) {
223 G4DynamicParticle* p = result->GetSecondary(0)->GetParticle();
224
225 if(p->GetKineticEnergy() > tcut) {
228 // G4cout << "recoil " << pdir << G4endl;
229 pdir.rotateUz(indir);
230 // G4cout << "recoil rotated " << pdir << G4endl;
231 p->SetMomentumDirection(pdir);
232
233 // in elastic scattering time and weight are not changed
234 G4Track* t = new G4Track(p, track.GetGlobalTime(),
235 track.GetPosition());
236 t->SetWeight(weight);
237 t->SetTouchableHandle(track.GetTouchableHandle());
239
240 } else {
241 edep += p->GetKineticEnergy();
242 delete p;
243 }
244 }
247 result->Clear();
248
249 return theTotalResult;
250}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
G4TrackStatus
@ fAlive
@ fStopAndKill
@ fStopButAlive
double G4double
Definition: G4Types.hh:83
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
const G4Element * SampleZandA(const G4DynamicParticle *, const G4Material *, G4Nucleus &target)
void SetMomentumDirection(const G4ThreeVector &aDirection)
const G4ThreeVector & GetMomentumDirection() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4String & GetName() const
Definition: G4Element.hh:126
G4double GetEnergyChange() const
void SetTrafoToLab(const G4LorentzRotation &aT)
G4double GetLocalEnergyDeposit() const
const G4ThreeVector & GetMomentumChange() const
std::size_t GetNumberOfSecondaries() const
G4HadSecondary * GetSecondary(size_t i)
G4DynamicParticle * GetParticle()
void Report(std::ostream &aS) const
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
const G4String & GetModelName() const
void SetRecoilEnergyThreshold(G4double val)
void FillResult(G4HadFinalState *aR, const G4Track &aT)
G4Nucleus * GetTargetNucleusPointer()
G4HadFinalState * CheckResult(const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
G4ParticleChange * theTotalResult
G4HadronicInteraction * ChooseHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, const G4Material *aMaterial, const G4Element *anElement)
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
G4CrossSectionDataStore * GetCrossSectionDataStore()
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
const G4String & GetName() const
Definition: G4Material.hh:175
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:651
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void AddSecondary(G4Track *aSecondary)
void ProposeEnergy(G4double finalEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
virtual void Initialize(const G4Track &)
G4ProcessManager * GetProcessManager() const
const G4String & GetParticleName() const
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
std::size_t size() const
const std::vector< G4double > * GetEnergyCutsVector(std::size_t pcIdx) const
static G4ProductionCutsTable * GetProductionCutsTable()
void SetWeight(G4double aValue)
void SetTouchableHandle(const G4TouchableHandle &apValue)
virtual G4double ComputeRatio(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)=0
void ProposeTrackStatus(G4TrackStatus status)
void ProposeNonIonizingEnergyDeposit(G4double anEnergyPart)
void ProposeWeight(G4double finalWeight)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetNumberOfSecondaries(G4int totSecondaries)
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:424
G4int verboseLevel
Definition: G4VProcess.hh:356

◆ ProcessDescription()

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

Reimplemented from G4HadronicProcess.

Definition at line 54 of file G4HadronElasticProcess.cc.

55{
56 outFile << "G4HadronElasticProcess handles the elastic scattering of \n"
57 << "hadrons by invoking the following hadronic model(s) and \n"
58 << "hadronic cross section(s).\n";
59}

◆ SetDiffraction()

void G4HadronElasticProcess::SetDiffraction ( G4HadronicInteraction hi,
G4VCrossSectionRatio xsr 
)

Definition at line 263 of file G4HadronElasticProcess.cc.

265{
266 if(hi && xsr) {
267 fDiffraction = hi;
268 fDiffractionRatio = xsr;
269 }
270}

Referenced by G4HadronHElasticPhysics::ConstructProcess().

◆ SetLowestEnergy()

void G4HadronElasticProcess::SetLowestEnergy ( G4double  )
virtual

Definition at line 252 of file G4HadronElasticProcess.cc.

253{
254 PrintWarning("G4HadronElasticProcess::SetLowestEnergy(..) ");
255}

◆ SetLowestEnergyNeutron()

void G4HadronElasticProcess::SetLowestEnergyNeutron ( G4double  )
virtual

Definition at line 258 of file G4HadronElasticProcess.cc.

259{
260 PrintWarning("G4HadronElasticProcess::SetLowestEnergyNeutron(..) ");
261}

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