Geant4  10.02.p03
G4HadronElasticProcess Class Reference

#include <G4HadronElasticProcess.hh>

Inheritance diagram for G4HadronElasticProcess:
Collaboration diagram for G4HadronElasticProcess:

Public Member Functions

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

Private Member Functions

G4HadronElasticProcessoperator= (const G4HadronElasticProcess &right)
 
 G4HadronElasticProcess (const G4HadronElasticProcess &)
 

Private Attributes

G4double lowestEnergy
 
G4bool isInitialised
 
G4HadronicInteractionfDiffraction
 
G4VCrossSectionRatiofDiffractionRatio
 

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, G4Material *aMaterial, 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 &)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4HadronicProcess
G4HadProjectile thePro
 
G4ParticleChange * theTotalResult
 
G4int epReportLevel
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChange * pParticleChange
 
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
 

Detailed Description

Definition at line 50 of file G4HadronElasticProcess.hh.

Constructor & Destructor Documentation

◆ G4HadronElasticProcess() [1/2]

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

Definition at line 50 of file G4HadronElasticProcess.cc.

53 {
55  lowestEnergy = 1.*keV;
56 }
G4VCrossSectionRatio * fDiffractionRatio
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
G4HadronicInteraction * fDiffraction
G4HadronicProcess(const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
static const double keV
Definition: G4SIunits.hh:213
Here is the call graph for this function:

◆ ~G4HadronElasticProcess()

G4HadronElasticProcess::~G4HadronElasticProcess ( )
virtual

Definition at line 58 of file G4HadronElasticProcess.cc.

59 {}

◆ G4HadronElasticProcess() [2/2]

G4HadronElasticProcess::G4HadronElasticProcess ( const G4HadronElasticProcess )
private

Member Function Documentation

◆ operator=()

G4HadronElasticProcess& G4HadronElasticProcess::operator= ( const G4HadronElasticProcess right)
private

◆ PostStepDoIt()

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

! is not needed for models inheriting G4HadronElastic

Reimplemented from G4HadronicProcess.

Definition at line 71 of file G4HadronElasticProcess.cc.

73 {
74  theTotalResult->Clear();
75  theTotalResult->Initialize(track);
76  G4double weight = track.GetWeight();
77  theTotalResult->ProposeWeight(weight);
78 
79  // For elastic scattering, _any_ result is considered an interaction
81 
82  G4double kineticEnergy = track.GetKineticEnergy();
83  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
84  const G4ParticleDefinition* part = dynParticle->GetDefinition();
85 
86  // NOTE: Very low energy scatters were causing numerical (FPE) errors
87  // in earlier releases; these limits have not been changed since.
88  if (kineticEnergy <= lowestEnergy) return theTotalResult;
89 
90  G4Material* material = track.GetMaterial();
91  G4Nucleus* targNucleus = GetTargetNucleusPointer();
92 
93  // Select element
94  G4Element* elm = 0;
95  try
96  {
97  elm = GetCrossSectionDataStore()->SampleZandA(dynParticle, material,
98  *targNucleus);
99  }
100  catch(G4HadronicException & aR)
101  {
103  aR.Report(ed);
104  DumpState(track,"SampleZandA",ed);
105  ed << " PostStepDoIt failed on element selection" << G4endl;
106  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had003",
107  FatalException, ed);
108  }
109 
110  // Initialize the hadronic projectile from the track
111  // G4cout << "track " << track.GetDynamicParticle()->Get4Momentum()<<G4endl;
112  G4HadProjectile theProj(track);
113  G4HadronicInteraction* hadi = 0;
114  G4HadFinalState* result = 0;
115 
116  if(fDiffraction)
117  {
118  G4double ratio =
119  fDiffractionRatio->ComputeRatio(part, kineticEnergy,
120  targNucleus->GetZ_asInt(),
121  targNucleus->GetA_asInt());
122  // diffraction is chosen
123  if(ratio > 0.0 && G4UniformRand() < ratio)
124  {
125  try
126  {
127 // if(part->GetParticleName() == "pi-")
128 // G4cout<<part->GetParticleName()<<"; "<<kineticEnergy/CLHEP::GeV<<" GeV; r = "<<ratio<<G4endl;
129 
130  result = fDiffraction->ApplyYourself(theProj, *targNucleus);
131  }
132  catch(G4HadronicException aR)
133  {
135  aR.Report(ed);
136  ed << "Call for " << fDiffraction->GetModelName() << G4endl;
137  ed << "Target element "<< elm->GetName()<<" Z= "
138  << targNucleus->GetZ_asInt()
139  << " A= " << targNucleus->GetA_asInt() << G4endl;
140  DumpState(track,"ApplyYourself",ed);
141  ed << " ApplyYourself failed" << G4endl;
142  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
143  FatalException, ed);
144  }
145  // Check the result for catastrophic energy non-conservation
146  result = CheckResult(theProj, *targNucleus, result);
147 
148  result->SetTrafoToLab(theProj.GetTrafoToLab());
150 
151  FillResult(result, track);
152 
153  if (epReportLevel != 0) {
154  CheckEnergyMomentumConservation(track, *targNucleus);
155  }
156  return theTotalResult;
157  }
158  }
159 
160  // ordinary elastic scattering
161  try
162  {
163  hadi = ChooseHadronicInteraction( theProj, *targNucleus, material, elm );
164  }
165  catch(G4HadronicException & aE)
166  {
168  aE.Report(ed);
169  ed << "Target element "<< elm->GetName()<<" Z= "
170  << targNucleus->GetZ_asInt() << " A= "
171  << targNucleus->GetA_asInt() << G4endl;
172  DumpState(track,"ChooseHadronicInteraction",ed);
173  ed << " No HadronicInteraction found out" << G4endl;
174  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had005",
175  FatalException, ed);
176  }
177 
178  size_t idx = track.GetMaterialCutsCouple()->GetIndex();
180  ->GetEnergyCutsVector(3)))[idx];
181  hadi->SetRecoilEnergyThreshold(tcut);
182 
183  if(verboseLevel>1) {
184  G4cout << "G4HadronElasticProcess::PostStepDoIt for "
185  << part->GetParticleName()
186  << " in " << material->GetName()
187  << " Target Z= " << targNucleus->GetZ_asInt()
188  << " A= " << targNucleus->GetA_asInt() << G4endl;
189  }
190 
191  try
192  {
193  result = hadi->ApplyYourself( theProj, *targNucleus);
194  }
195  catch(G4HadronicException aR)
196  {
198  aR.Report(ed);
199  ed << "Call for " << hadi->GetModelName() << G4endl;
200  ed << "Target element "<< elm->GetName()<<" Z= "
201  << targNucleus->GetZ_asInt()
202  << " A= " << targNucleus->GetA_asInt() << G4endl;
203  DumpState(track,"ApplyYourself",ed);
204  ed << " ApplyYourself failed" << G4endl;
205  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
206  FatalException, ed);
207  }
208 
209  // Check the result for catastrophic energy non-conservation
210  // cannot be applied because is not guranteed that recoil
211  // nucleus is created
212  // result = CheckResult(theProj, targNucleus, result);
213 
214  // directions
215  G4ThreeVector indir = track.GetMomentumDirection();
217  G4ThreeVector it(0., 0., 1.);
218  G4ThreeVector outdir = result->GetMomentumChange();
219 
220  if(verboseLevel>1) {
221  G4cout << "Efin= " << result->GetEnergyChange()
222  << " de= " << result->GetLocalEnergyDeposit()
223  << " nsec= " << result->GetNumberOfSecondaries()
224  << " dir= " << outdir
225  << G4endl;
226  }
227 
228  // energies
229  G4double edep = result->GetLocalEnergyDeposit();
230  G4double efinal = result->GetEnergyChange();
231  if(efinal < 0.0) { efinal = 0.0; }
232  if(edep < 0.0) { edep = 0.0; }
233 
234  // NOTE: Very low energy scatters were causing numerical (FPE) errors
235  // in earlier releases; these limits have not been changed since.
236  if(efinal <= lowestEnergy) {
237  edep += efinal;
238  efinal = 0.0;
239  }
240 
241  // primary change
242  theTotalResult->ProposeEnergy(efinal);
243 
244  G4TrackStatus status = track.GetTrackStatus();
245  if(efinal > 0.0) {
246  outdir.rotate(phi, it);
247  outdir.rotateUz(indir);
248  theTotalResult->ProposeMomentumDirection(outdir);
249  } else {
250  if(part->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
251  { status = fStopButAlive; }
252  else { status = fStopAndKill; }
253  theTotalResult->ProposeTrackStatus(status);
254  }
255 
256  //G4cout << "Efinal= " << efinal << " TrackStatus= " << status << G4endl;
257 
258  theTotalResult->SetNumberOfSecondaries(0);
259 
260  // recoil
261  if(result->GetNumberOfSecondaries() > 0) {
262  G4DynamicParticle* p = result->GetSecondary(0)->GetParticle();
263 
264  if(p->GetKineticEnergy() > tcut) {
265  theTotalResult->SetNumberOfSecondaries(1);
267  // G4cout << "recoil " << pdir << G4endl;
269  pdir.rotate(phi, it);
270  pdir.rotateUz(indir);
271  // G4cout << "recoil rotated " << pdir << G4endl;
272  p->SetMomentumDirection(pdir);
273 
274  // in elastic scattering time and weight are not changed
275  G4Track* t = new G4Track(p, track.GetGlobalTime(),
276  track.GetPosition());
277  t->SetWeight(weight);
278  t->SetTouchableHandle(track.GetTouchableHandle());
279  theTotalResult->AddSecondary(t);
280 
281  } else {
282  edep += p->GetKineticEnergy();
283  delete p;
284  }
285  }
286  theTotalResult->ProposeLocalEnergyDeposit(edep);
287  theTotalResult->ProposeNonIonizingEnergyDeposit(edep);
288  result->Clear();
289 
290  return theTotalResult;
291 }
virtual G4double ComputeRatio(const G4ParticleDefinition *, G4double kinEnergy, G4int, G4int)=0
G4int verboseLevel
Definition: G4VProcess.hh:368
G4int GetNumberOfSecondaries() const
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4HadSecondary * GetSecondary(size_t i)
G4double GetEnergyChange() const
void SetMomentumDirection(const G4ThreeVector &aDirection)
G4ProcessManager * GetProcessManager() const
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
double weight
Definition: plottest35.C:25
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
Double_t edep
G4VCrossSectionRatio * fDiffractionRatio
string material
Definition: eplot.py:19
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
void SetTrafoToLab(const G4LorentzRotation &aT)
G4double GetKineticEnergy() const
void FillResult(G4HadFinalState *aR, const G4Track &aT)
const G4String & GetParticleName() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4double GetLocalEnergyDeposit() const
G4GLOB_DLL std::ostream G4cout
G4ParticleChange * theTotalResult
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:38
G4CrossSectionDataStore * GetCrossSectionDataStore()
TString part[npart]
const G4ThreeVector & GetMomentumChange() const
G4int size() const
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
const G4String & GetName() const
Definition: G4Element.hh:127
static G4ProductionCutsTable * GetProductionCutsTable()
G4Nucleus * GetTargetNucleusPointer()
const G4ThreeVector & GetMomentumDirection() const
G4DynamicParticle * GetParticle()
G4HadronicInteraction * fDiffraction
G4HadronicInteraction * ChooseHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, G4Material *aMaterial, G4Element *anElement)
G4ParticleDefinition * GetDefinition() const
#define G4endl
Definition: G4ios.hh:61
void SetRecoilEnergyThreshold(G4double val)
G4Element * SampleZandA(const G4DynamicParticle *, const G4Material *, G4Nucleus &target)
double G4double
Definition: G4Types.hh:76
const G4String & GetModelName() const
static const double twopi
Definition: SystemOfUnits.h:54
const G4String & GetName() const
Definition: G4Material.hh:178
Hep3Vector & rotate(double, const Hep3Vector &)
Definition: ThreeVectorR.cc:28
void Report(std::ostream &aS)
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4HadFinalState * CheckResult(const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
Here is the call graph for this function:

◆ PreparePhysicsTable()

void G4HadronElasticProcess::PreparePhysicsTable ( const G4ParticleDefinition part)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 294 of file G4HadronElasticProcess.cc.

295 {
296  if(!isInitialised) {
297  isInitialised = true;
298  if(G4Neutron::Neutron() == &part) { lowestEnergy = 1.e-6*eV; }
299  }
301 }
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static const double eV
Definition: G4SIunits.hh:212
Here is the call graph for this function:

◆ ProcessDescription()

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

Reimplemented from G4HadronicProcess.

Definition at line 61 of file G4HadronElasticProcess.cc.

62 {
63 
64  outFile << "G4HadronElasticProcess handles the elastic scattering of \n"
65  << "hadrons by invoking the following hadronic model(s) and \n"
66  << "hadronic cross section(s).\n";
67 
68 }

◆ SetDiffraction()

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

Definition at line 315 of file G4HadronElasticProcess.cc.

317 {
318  if(hi && xsr) {
319  fDiffraction = hi;
320  fDiffractionRatio = xsr;
321  }
322 }
G4VCrossSectionRatio * fDiffractionRatio
G4HadronicInteraction * fDiffraction

◆ SetLowestEnergy()

void G4HadronElasticProcess::SetLowestEnergy ( G4double  val)
virtual

Definition at line 303 of file G4HadronElasticProcess.cc.

304 {
305  lowestEnergy = val;
306 }

◆ SetLowestEnergyNeutron()

void G4HadronElasticProcess::SetLowestEnergyNeutron ( G4double  val)
virtual

Definition at line 309 of file G4HadronElasticProcess.cc.

310 {
311  lowestEnergy = val;
312  G4HadronicDeprecate("G4HadronElasticProcess::SetLowestEnergyNeutron()");
313 }
#define G4HadronicDeprecate(name)

Member Data Documentation

◆ fDiffraction

G4HadronicInteraction* G4HadronElasticProcess::fDiffraction
private

Definition at line 83 of file G4HadronElasticProcess.hh.

◆ fDiffractionRatio

G4VCrossSectionRatio* G4HadronElasticProcess::fDiffractionRatio
private

Definition at line 84 of file G4HadronElasticProcess.hh.

◆ isInitialised

G4bool G4HadronElasticProcess::isInitialised
private

Definition at line 82 of file G4HadronElasticProcess.hh.

◆ lowestEnergy

G4double G4HadronElasticProcess::lowestEnergy
private

Definition at line 81 of file G4HadronElasticProcess.hh.


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