Geant4 11.1.1
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)
 
void StartTracking (G4Track *track) override
 
G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double, G4ForceCondition *) override
 
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 &)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *) override
 
const G4NucleusGetTargetNucleus () const
 
G4NucleusGetTargetNucleusPointer ()
 
const G4IsotopeGetTargetIsotope ()
 
G4double ComputeCrossSection (const G4ParticleDefinition *, const G4Material *, const G4double kinEnergy)
 
G4HadXSType CrossSectionType () const
 
void SetCrossSectionType (G4HadXSType val)
 
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 ()
 
std::vector< G4TwoPeaksHadXS * > * TwoPeaksXS () const
 
std::vector< G4double > * EnergyOfCrossSectionMax () const
 
G4HadronicProcessoperator= (const G4HadronicProcess &right)=delete
 
 G4HadronicProcess (const G4HadronicProcess &)=delete
 
- 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 &)
 
virtual G4double GetCrossSection (const G4double, const G4MaterialCutsCouple *)
 
virtual G4double MinPrimaryEnergy (const G4ParticleDefinition *, const G4Material *)
 
- 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 &)
 

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)
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
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
 
G4CrossSectionDataStoretheCrossSectionDataStore
 
G4double fWeight = 1.0
 
G4double aScaleFactor = 1.0
 
G4double theLastCrossSection = 0.0
 
G4double mfpKinEnergy = DBL_MAX
 
G4long epReportLevel = 0
 
G4HadXSType fXSType = fHadNoIntegral
 
- 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 const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
74 G4double kineticEnergy = dynParticle->GetKineticEnergy();
75 G4TrackStatus status = track.GetTrackStatus();
76 if(kineticEnergy == 0.0 || track.GetTrackStatus() != fAlive) {
77 return theTotalResult;
78 }
79
80 const G4Material* material = track.GetMaterial();
81
82 // check only for charged particles
83 if(fXSType != fHadNoIntegral) {
86 theCrossSectionDataStore->ComputeCrossSection(dynParticle, material);
88 // No interaction
89 return theTotalResult;
90 }
91 }
92
93 const G4ParticleDefinition* part = dynParticle->GetDefinition();
94 G4Nucleus* targNucleus = GetTargetNucleusPointer();
95
96 // Select element
97 const G4Element* elm =
98 theCrossSectionDataStore->SampleZandA(dynParticle, material, *targNucleus);
99
100 // Initialize the hadronic projectile from the track
101 G4HadProjectile theProj(track);
102 G4HadronicInteraction* hadi = nullptr;
103 G4HadFinalState* result = nullptr;
104
105 if(nullptr != fDiffraction) {
106 G4double ratio =
107 fDiffractionRatio->ComputeRatio(part, kineticEnergy,
108 targNucleus->GetZ_asInt(),
109 targNucleus->GetA_asInt());
110 // diffraction is chosen
111 if(ratio > 0.0 && G4UniformRand() < ratio)
112 {
113 try
114 {
115 result = fDiffraction->ApplyYourself(theProj, *targNucleus);
116 }
117 catch(G4HadronicException & aR)
118 {
120 aR.Report(ed);
121 ed << "Call for " << fDiffraction->GetModelName() << G4endl;
122 ed << part->GetParticleName()
123 << " off target element " << elm->GetName() << " Z= "
124 << targNucleus->GetZ_asInt()
125 << " A= " << targNucleus->GetA_asInt() << G4endl;
126 DumpState(track,"ApplyYourself",ed);
127 ed << " ApplyYourself failed" << G4endl;
128 G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
129 FatalException, ed);
130 }
131 // Check the result for catastrophic energy non-conservation
132 result = CheckResult(theProj, *targNucleus, result);
133 result->SetTrafoToLab(theProj.GetTrafoToLab());
134
135 // The following method of the base class takes care also of setting
136 // the creator model ID for the secondaries that are created
137 FillResult(result, track);
138
139 if (epReportLevel != 0) {
140 CheckEnergyMomentumConservation(track, *targNucleus);
141 }
142 return theTotalResult;
143 }
144 }
145
146 // ordinary elastic scattering
147 hadi = ChooseHadronicInteraction( theProj, *targNucleus, material, elm );
148 if(nullptr == hadi) {
150 ed << part->GetParticleName()
151 << " off target element " << elm->GetName() << " Z= "
152 << targNucleus->GetZ_asInt() << " A= "
153 << targNucleus->GetA_asInt() << G4endl;
154 DumpState(track,"ChooseHadronicInteraction",ed);
155 ed << " No HadronicInteraction found out" << G4endl;
156 G4Exception("G4HadronElasticProcess::PostStepDoIt", "had005",
157 FatalException, ed);
158 return theTotalResult;
159 }
160
161 size_t idx = track.GetMaterialCutsCouple()->GetIndex();
163 ->GetEnergyCutsVector(3)))[idx];
164 hadi->SetRecoilEnergyThreshold(tcut);
165 /*
166 if(verboseLevel>1) {
167 G4cout << "G4HadronElasticProcess::PostStepDoIt for "
168 << part->GetParticleName()
169 << " in " << material->GetName()
170 << " Target Z= " << targNucleus->GetZ_asInt()
171 << " A= " << targNucleus->GetA_asInt()
172 << " Tcut(MeV)= " << tcut << G4endl;
173 }
174 */
175 result = hadi->ApplyYourself( theProj, *targNucleus);
176
177 // Check the result for catastrophic energy non-conservation
178 // cannot be applied because is not guranteed that recoil
179 // nucleus is created
180 // result = CheckResult(theProj, targNucleus, result);
181
182 // directions
183 G4ThreeVector indir = track.GetMomentumDirection();
184 G4ThreeVector outdir = result->GetMomentumChange();
185 /*
186 if(verboseLevel>1) {
187 G4cout << "Efin= " << result->GetEnergyChange()
188 << " de= " << result->GetLocalEnergyDeposit()
189 << " nsec= " << result->GetNumberOfSecondaries()
190 << " dir= " << outdir
191 << G4endl;
192 }
193 */
194 // energies
195 G4double edep = std::max(result->GetLocalEnergyDeposit(), 0.0);
196 G4double efinal = std::max(result->GetEnergyChange(), 0.0);
197
198 // primary change
200
201 if(efinal > 0.0) {
202 outdir.rotateUz(indir);
204 } else {
205 if(part->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
206 { status = fStopButAlive; }
207 else { status = fStopAndKill; }
209 }
210 /*
211 G4cout << "Efinal= " << efinal << " TrackStatus= " << status
212 << " time(ns)=" << track.GetGlobalTime()/ns << G4endl;
213 */
215
216 // recoil
217 if(result->GetNumberOfSecondaries() > 0) {
218 G4DynamicParticle* p = result->GetSecondary(0)->GetParticle();
219
220 if(p->GetKineticEnergy() > tcut) {
223 // G4cout << "recoil " << pdir << G4endl;
224 pdir.rotateUz(indir);
225 // G4cout << "recoil rotated " << pdir << G4endl;
226 p->SetMomentumDirection(pdir);
227
228 // in elastic scattering time and weight are not changed
229 G4Track* t = new G4Track(p, track.GetGlobalTime(),
230 track.GetPosition());
231 t->SetWeight(weight);
232 t->SetTouchableHandle(track.GetTouchableHandle());
233 G4int secID = G4PhysicsModelCatalog::GetModelID( "model_" + hadi->GetModelName() );
234 if ( secID > 0 ) t->SetCreatorModelID(secID);
236
237 } else {
238 edep += p->GetKineticEnergy();
239 delete p;
240 }
241 }
244 result->Clear();
245
246 return theTotalResult;
247}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
@ fHadNoIntegral
Definition: G4HadXSTypes.hh:45
G4TrackStatus
@ fAlive
@ fStopAndKill
@ fStopButAlive
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
G4double ComputeCrossSection(const G4DynamicParticle *, const G4Material *)
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:127
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)
G4CrossSectionDataStore * theCrossSectionDataStore
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:691
G4int GetA_asInt() const
Definition: G4Nucleus.hh:99
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:105
void AddSecondary(G4Track *aSecondary)
void Initialize(const G4Track &) override
void ProposeEnergy(G4double finalEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4ProcessManager * GetProcessManager() const
const G4String & GetParticleName() const
static G4int GetModelID(const G4int modelIndex)
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)
void SetCreatorModelID(const G4int id)
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)
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:335
#define DBL_MAX
Definition: templates.hh:62

◆ 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 260 of file G4HadronElasticProcess.cc.

262{
263 if(hi && xsr) {
264 fDiffraction = hi;
265 fDiffractionRatio = xsr;
266 }
267}

Referenced by G4HadronHElasticPhysics::ConstructProcess().

◆ SetLowestEnergy()

void G4HadronElasticProcess::SetLowestEnergy ( G4double  )
virtual

Definition at line 249 of file G4HadronElasticProcess.cc.

250{
251 PrintWarning("G4HadronElasticProcess::SetLowestEnergy(..) ");
252}

◆ SetLowestEnergyNeutron()

void G4HadronElasticProcess::SetLowestEnergyNeutron ( G4double  )
virtual

Definition at line 255 of file G4HadronElasticProcess.cc.

256{
257 PrintWarning("G4HadronElasticProcess::SetLowestEnergyNeutron(..) ");
258}

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