Geant4.10
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions
G4HadronElasticProcess Class Reference

#include <G4HadronElasticProcess.hh>

Inheritance diagram for G4HadronElasticProcess:
G4HadronicProcess G4VDiscreteProcess G4VProcess

Public Member Functions

 G4HadronElasticProcess (const G4String &procName="hadElastic")
 
virtual ~G4HadronElasticProcess ()
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual void SetLowestEnergy (G4double)
 
virtual void SetLowestEnergyNeutron (G4double)
 
virtual void Description () const
 
- 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)
 
G4EnergyRangeManagerGetManagerPointer ()
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
 
const G4NucleusGetTargetNucleus () const
 
const G4IsotopeGetTargetIsotope ()
 
virtual void ProcessDescription (std::ostream &outFile) const
 
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 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 ()
 
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 &)
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Member Functions inherited from G4HadronicProcess
G4HadronicInteractionChooseHadronicInteraction (G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement)
 
G4NucleusGetTargetNucleusPointer ()
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
const G4EnergyRangeManagerGetEnergyRangeManager () const
 
void SetEnergyRangeManager (const G4EnergyRangeManager &value)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
G4HadFinalStateCheckResult (const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result) const
 
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
 
G4ParticleChangetheTotalResult
 
G4int epReportLevel
 
- 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
 

Detailed Description

Definition at line 49 of file G4HadronElasticProcess.hh.

Constructor & Destructor Documentation

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

Definition at line 48 of file G4HadronElasticProcess.cc.

References G4HadronicProcess::AddDataSet(), and python.hepunit::keV.

49  : G4HadronicProcess(pName, fHadronElastic), isInitialised(false)
50 {
52  lowestEnergy = 1.*keV;
53 }
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
G4HadronicProcess(const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
G4HadronElasticProcess::~G4HadronElasticProcess ( )
virtual

Definition at line 55 of file G4HadronElasticProcess.cc.

56 {}

Member Function Documentation

void G4HadronElasticProcess::Description ( ) const
virtual

Definition at line 58 of file G4HadronElasticProcess.cc.

References G4VProcess::GetProcessName(), and outFile.

59 {
60  char* dirName = getenv("G4PhysListDocDir");
61  if (dirName) {
62  std::ofstream outFile;
63  G4String outFileName = GetProcessName() + ".html";
64  G4String pathName = G4String(dirName) + "/" + outFileName;
65  outFile.open(pathName);
66  outFile << "<html>\n";
67  outFile << "<head>\n";
68 
69  outFile << "<title>Description of G4HadronElasticProcess</title>\n";
70  outFile << "</head>\n";
71  outFile << "<body>\n";
72 
73  outFile << "G4HadronElasticProcess handles the elastic scattering of\n"
74  << "hadrons by invoking one or more hadronic models and one or\n"
75  << "more hadronic cross sections.\n";
76 
77  outFile << "</body>\n";
78  outFile << "</html>\n";
79  outFile.close();
80  }
81 }
std::ofstream outFile
Definition: GammaRayTel.cc:68
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
G4VParticleChange * G4HadronElasticProcess::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
virtual

! is not needed for models inheriting G4HadronElastic

Reimplemented from G4HadronicProcess.

Definition at line 84 of file G4HadronElasticProcess.cc.

References G4ParticleChange::AddSecondary(), G4HadronicInteraction::ApplyYourself(), G4HadronicProcess::ChooseHadronicInteraction(), G4HadFinalState::Clear(), G4VParticleChange::Clear(), G4VProcess::ClearNumberOfInteractionLengthLeft(), G4HadronicProcess::DumpState(), FatalException, fStopAndKill, fStopButAlive, G4cout, G4endl, G4Exception(), G4UniformRand, G4Nucleus::GetA_asInt(), G4ProcessManager::GetAtRestProcessVector(), G4HadronicProcess::GetCrossSectionDataStore(), G4DynamicParticle::GetDefinition(), G4Track::GetDynamicParticle(), G4HadFinalState::GetEnergyChange(), G4ProductionCutsTable::GetEnergyCutsVector(), G4Track::GetGlobalTime(), G4MaterialCutsCouple::GetIndex(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetKineticEnergy(), G4HadFinalState::GetLocalEnergyDeposit(), G4Track::GetMaterial(), G4Track::GetMaterialCutsCouple(), G4HadronicInteraction::GetModelName(), G4HadFinalState::GetMomentumChange(), G4DynamicParticle::GetMomentumDirection(), G4Track::GetMomentumDirection(), G4Element::GetName(), G4Material::GetName(), G4HadFinalState::GetNumberOfSecondaries(), G4HadSecondary::GetParticle(), G4ParticleDefinition::GetParticleName(), G4Track::GetPosition(), G4ParticleDefinition::GetProcessManager(), G4ProductionCutsTable::GetProductionCutsTable(), G4HadFinalState::GetSecondary(), G4HadronicProcess::GetTargetNucleusPointer(), G4Track::GetTouchableHandle(), G4Track::GetTrackStatus(), G4Track::GetWeight(), G4Nucleus::GetZ_asInt(), G4ParticleChange::Initialize(), eplot::material, G4ParticleChange::ProposeEnergy(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChange::ProposeMomentumDirection(), G4VParticleChange::ProposeNonIonizingEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), G4VParticleChange::ProposeWeight(), G4HadronicException::Report(), CLHEP::Hep3Vector::rotate(), CLHEP::Hep3Vector::rotateUz(), G4CrossSectionDataStore::SampleZandA(), G4DynamicParticle::SetMomentumDirection(), G4VParticleChange::SetNumberOfSecondaries(), G4HadronicInteraction::SetRecoilEnergyThreshold(), G4Track::SetTouchableHandle(), G4Track::SetWeight(), G4ProcessVector::size(), status, G4HadronicProcess::theTotalResult, and G4VProcess::verboseLevel.

86 {
88  theTotalResult->Initialize(track);
89  G4double weight = track.GetWeight();
91 
92  // For elastic scattering, _any_ result is considered an interaction
94 
95  G4double kineticEnergy = track.GetKineticEnergy();
96  const G4DynamicParticle* dynParticle = track.GetDynamicParticle();
97  const G4ParticleDefinition* part = dynParticle->GetDefinition();
98 
99  // NOTE: Very low energy scatters were causing numerical (FPE) errors
100  // in earlier releases; these limits have not been changed since.
101  if (kineticEnergy <= lowestEnergy) return theTotalResult;
102 
103  G4Material* material = track.GetMaterial();
104  G4Nucleus* targNucleus = GetTargetNucleusPointer();
105 
106  // Select element
107  G4Element* elm = 0;
108  try
109  {
110  elm = GetCrossSectionDataStore()->SampleZandA(dynParticle, material,
111  *targNucleus);
112  }
113  catch(G4HadronicException & aR)
114  {
116  aR.Report(ed);
117  DumpState(track,"SampleZandA",ed);
118  ed << " PostStepDoIt failed on element selection" << G4endl;
119  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had003",
120  FatalException, ed);
121  }
122  G4HadronicInteraction* hadi = 0;
123  try
124  {
125  hadi = ChooseHadronicInteraction( kineticEnergy, material, elm );
126  }
127  catch(G4HadronicException & aE)
128  {
130  aE.Report(ed);
131  ed << "Target element "<< elm->GetName()<<" Z= "
132  << targNucleus->GetZ_asInt() << " A= "
133  << targNucleus->GetA_asInt() << G4endl;
134  DumpState(track,"ChooseHadronicInteraction",ed);
135  ed << " No HadronicInteraction found out" << G4endl;
136  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had005",
137  FatalException, ed);
138  }
139 
140  size_t idx = track.GetMaterialCutsCouple()->GetIndex();
142  ->GetEnergyCutsVector(3)))[idx];
143  hadi->SetRecoilEnergyThreshold(tcut);
144 
145  // Initialize the hadronic projectile from the track
146  // G4cout << "track " << track.GetDynamicParticle()->Get4Momentum()<<G4endl;
147  G4HadProjectile theProj(track);
148  if(verboseLevel>1) {
149  G4cout << "G4HadronElasticProcess::PostStepDoIt for "
150  << part->GetParticleName()
151  << " in " << material->GetName()
152  << " Target Z= " << targNucleus->GetZ_asInt()
153  << " A= " << targNucleus->GetA_asInt() << G4endl;
154  }
155 
156  G4HadFinalState* result = 0;
157  try
158  {
159  result = hadi->ApplyYourself( theProj, *targNucleus);
160  }
161  catch(G4HadronicException aR)
162  {
164  aR.Report(ed);
165  ed << "Call for " << hadi->GetModelName() << G4endl;
166  ed << "Target element "<< elm->GetName()<<" Z= "
167  << targNucleus->GetZ_asInt()
168  << " A= " << targNucleus->GetA_asInt() << G4endl;
169  DumpState(track,"ApplyYourself",ed);
170  ed << " ApplyYourself failed" << G4endl;
171  G4Exception("G4HadronElasticProcess::PostStepDoIt", "had006",
172  FatalException, ed);
173  }
174 
175  // Check the result for catastrophic energy non-conservation
176  // cannot be applied because is not guranteed that recoil
177  // nucleus is created
178  // result = CheckResult(theProj, targNucleus, result);
179 
180  // directions
181  G4ThreeVector indir = track.GetMomentumDirection();
182  G4double phi = CLHEP::twopi*G4UniformRand();
183  G4ThreeVector it(0., 0., 1.);
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 = result->GetLocalEnergyDeposit();
196  G4double efinal = result->GetEnergyChange();
197  if(efinal < 0.0) { efinal = 0.0; }
198  if(edep < 0.0) { edep = 0.0; }
199 
200  // NOTE: Very low energy scatters were causing numerical (FPE) errors
201  // in earlier releases; these limits have not been changed since.
202  if(efinal <= lowestEnergy) {
203  edep += efinal;
204  efinal = 0.0;
205  }
206 
207  // primary change
208  theTotalResult->ProposeEnergy(efinal);
209 
210  G4TrackStatus status = track.GetTrackStatus();
211  if(efinal > 0.0) {
212  outdir.rotate(phi, it);
213  outdir.rotateUz(indir);
215  } else {
216  if(part->GetProcessManager()->GetAtRestProcessVector()->size() > 0)
217  { status = fStopButAlive; }
218  else { status = fStopAndKill; }
220  }
221 
222  //G4cout << "Efinal= " << efinal << " TrackStatus= " << status << G4endl;
223 
225 
226  // recoil
227  if(result->GetNumberOfSecondaries() > 0) {
228  G4DynamicParticle* p = result->GetSecondary(0)->GetParticle();
229 
230  if(p->GetKineticEnergy() > tcut) {
233  // G4cout << "recoil " << pdir << G4endl;
234  //!! is not needed for models inheriting G4HadronElastic
235  pdir.rotate(phi, it);
236  pdir.rotateUz(indir);
237  // G4cout << "recoil rotated " << pdir << G4endl;
238  p->SetMomentumDirection(pdir);
239 
240  // in elastic scattering time and weight are not changed
241  G4Track* t = new G4Track(p, track.GetGlobalTime(),
242  track.GetPosition());
243  t->SetWeight(weight);
244  t->SetTouchableHandle(track.GetTouchableHandle());
246 
247  } else {
248  edep += p->GetKineticEnergy();
249  delete p;
250  }
251  }
254  result->Clear();
255 
256  return theTotalResult;
257 }
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
const std::vector< G4double > * GetEnergyCutsVector(size_t pcIdx) const
G4int verboseLevel
Definition: G4VProcess.hh:368
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4HadSecondary * GetSecondary(size_t i)
G4double GetKineticEnergy() const
const char * p
Definition: xmltok.h:285
const G4String & GetName() const
Definition: G4Material.hh:176
void SetTouchableHandle(const G4TouchableHandle &apValue)
void SetMomentumDirection(const G4ThreeVector &aDirection)
G4double GetEnergyChange() const
G4ParticleDefinition * GetDefinition() const
const G4String & GetModelName() const
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
G4ProcessManager * GetProcessManager() const
const G4String & GetParticleName() const
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
string material
Definition: eplot.py:19
void ProposeWeight(G4double finalWeight)
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
G4ProcessVector * GetAtRestProcessVector(G4ProcessVectorTypeIndex typ=typeGPIL) const
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
G4ParticleChange * theTotalResult
void ProposeNonIonizingEnergyDeposit(G4double anEnergyPart)
const G4ThreeVector & GetMomentumDirection() const
G4HadronicInteraction * ChooseHadronicInteraction(G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:72
G4CrossSectionDataStore * GetCrossSectionDataStore()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
int status
Definition: tracer.cxx:24
G4int size() const
virtual void Initialize(const G4Track &)
static G4ProductionCutsTable * GetProductionCutsTable()
G4Nucleus * GetTargetNucleusPointer()
const G4ThreeVector & GetMomentumChange() const
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * GetParticle()
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void ProposeEnergy(G4double finalEnergy)
void AddSecondary(G4Track *aSecondary)
#define G4endl
Definition: G4ios.hh:61
void SetRecoilEnergyThreshold(G4double val)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
G4Element * SampleZandA(const G4DynamicParticle *, const G4Material *, G4Nucleus &target)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
const G4String & GetName() const
Definition: G4Element.hh:127
G4TrackStatus
Hep3Vector & rotate(double, const Hep3Vector &)
Definition: ThreeVectorR.cc:28
void Report(std::ostream &aS)
G4int GetNumberOfSecondaries() const
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4double GetLocalEnergyDeposit() const
void G4HadronElasticProcess::PreparePhysicsTable ( const G4ParticleDefinition part)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 260 of file G4HadronElasticProcess.cc.

References python.hepunit::eV, G4Neutron::Neutron(), and G4HadronicProcess::PreparePhysicsTable().

261 {
262  if(!isInitialised) {
263  isInitialised = true;
264  if(G4Neutron::Neutron() == &part) { lowestEnergy = 1.e-6*eV; }
265  }
267 }
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
void G4HadronElasticProcess::SetLowestEnergy ( G4double  val)
virtual

Definition at line 269 of file G4HadronElasticProcess.cc.

270 {
271  lowestEnergy = val;
272 }
void G4HadronElasticProcess::SetLowestEnergyNeutron ( G4double  val)
virtual

Definition at line 275 of file G4HadronElasticProcess.cc.

References G4HadronicDeprecate.

276 {
277  lowestEnergy = val;
278  G4HadronicDeprecate("G4HadronElasticProcess::SetLowestEnergyNeutron()");
279 }
#define G4HadronicDeprecate(name)

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