Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4HadronStoppingProcess Class Reference

#include <G4HadronStoppingProcess.hh>

Inheritance diagram for G4HadronStoppingProcess:
Collaboration diagram for G4HadronStoppingProcess:

Public Member Functions

 G4HadronStoppingProcess (const G4String &name="hadronCaptureAtRest")
 
virtual ~G4HadronStoppingProcess ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual void ProcessDescription (std::ostream &outFile) const
 
void SetElementSelector (G4ElementSelector *ptr)
 
void SetEmCascade (G4HadronicInteraction *ptr)
 
void SetBoundDecay (G4HadronicInteraction *ptr)
 
- 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=nullptr)
 
G4double GetMicroscopicCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
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 SetIntegral (G4bool val)
 
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 AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
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 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 &)
 

Protected Member Functions

G4double GetMeanLifeTime (const G4Track &, G4ForceCondition *)
 
- 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 ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- 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 63 of file G4HadronStoppingProcess.hh.

Constructor & Destructor Documentation

G4HadronStoppingProcess::G4HadronStoppingProcess ( const G4String name = "hadronCaptureAtRest")
explicit

Definition at line 64 of file G4HadronStoppingProcess.cc.

66  fElementSelector(new G4ElementSelector()),
67  fEmCascade(new G4EmCaptureCascade()), // Owned by InteractionRegistry
68  fBoundDecay(0),
69  emcID(-1),
70  ncID(-1),
71  dioID(-1)
72 {
73  // Modify G4VProcess flags to emulate G4VRest instead of G4VDiscrete
74  enableAtRestDoIt = true;
75  enablePostStepDoIt = false;
76 
78 }
static G4HadronicProcessStore * Instance()
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
void RegisterExtraProcess(G4VProcess *)
G4HadronicProcess(const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)

Here is the call graph for this function:

G4HadronStoppingProcess::~G4HadronStoppingProcess ( )
virtual

Definition at line 82 of file G4HadronStoppingProcess.cc.

83 {
84  //G4HadronicProcessStore::Instance()->DeRegisterExtraProcess(this);
85  delete fElementSelector;
86  // NOTE: fEmCascade and fEmBoundDecay owned by registry, not locally
87 }

Member Function Documentation

G4VParticleChange * G4HadronStoppingProcess::AtRestDoIt ( const G4Track track,
const G4Step  
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 134 of file G4HadronStoppingProcess.cc.

136 {
137  // if primary is not Alive then do nothing
138  theTotalResult->Initialize(track);
139 
140  G4Nucleus* nucleus = GetTargetNucleusPointer();
141  G4Element* elm = fElementSelector->SelectZandA(track, nucleus);
142 
143  G4HadFinalState* result = 0;
144  thePro.Initialise(track);
145 
146  // save track time an dstart capture from zero time
147  thePro.SetGlobalTime(0.0);
148  G4double time0 = track.GetGlobalTime();
149 
150  G4bool nuclearCapture = true;
151 
152  // Do the electromagnetic cascade in the nuclear field.
153  // EM cascade should keep G4HadFinalState object,
154  // because it will not be deleted at the end of this method
155  //
156  result = fEmCascade->ApplyYourself(thePro, *nucleus);
157  G4double ebound = result->GetLocalEnergyDeposit();
158  G4double edep = 0.0;
159  G4int nSecondaries = result->GetNumberOfSecondaries();
160  G4int nEmCascadeSec = nSecondaries;
161 
162  // Try decay from bound level
163  // For mu- the time of projectile should be changed.
164  // Decay should keep G4HadFinalState object,
165  // because it will not be deleted at the end of this method.
166  //
167  thePro.SetBoundEnergy(ebound);
168  if(fBoundDecay) {
169  G4HadFinalState* resultDecay =
170  fBoundDecay->ApplyYourself(thePro, *nucleus);
171  G4int n = resultDecay->GetNumberOfSecondaries();
172  if(0 < n) {
173  nSecondaries += n;
174  result->AddSecondaries(resultDecay);
175  }
176  if(resultDecay->GetStatusChange() == stopAndKill) {
177  nuclearCapture = false;
178  }
179  resultDecay->Clear();
180  }
181 
182  if(nuclearCapture) {
183 
184  // delay of capture
185  G4double capTime = thePro.GetGlobalTime();
186  thePro.SetGlobalTime(0.0);
187 
188  // select model
190  try {
191  model = ChooseHadronicInteraction(thePro, *nucleus,
192  track.GetMaterial(), elm);
193  }
194  catch(G4HadronicException & aE) {
196  ed << "Target element "<<elm->GetName()<<" Z= "
197  << nucleus->GetZ_asInt() << " A= "
198  << nucleus->GetA_asInt() << G4endl;
199  DumpState(track,"ChooseHadronicInteraction",ed);
200  ed << " No HadronicInteraction found out" << G4endl;
201  G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had005",
202  FatalException, ed);
203  }
204 
205  G4HadFinalState* resultNuc = 0;
206  G4int reentryCount = 0;
207  do {
208  // sample final state
209  // nuclear interaction should keep G4HadFinalState object
210  // model should define time of each secondary particle
211  try {
212  resultNuc = model->ApplyYourself(thePro, *nucleus);
213  ++reentryCount;
214  }
215  catch(G4HadronicException aR) {
217  ed << "Call for " << model->GetModelName() << G4endl;
218  ed << "Target element "<<elm->GetName()<<" Z= "
219  << nucleus->GetZ_asInt()
220  << " A= " << nucleus->GetA_asInt() << G4endl;
221  DumpState(track,"ApplyYourself",ed);
222  ed << " ApplyYourself failed" << G4endl;
223  G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006",
224  FatalException, ed);
225  }
226 
227  // Check the result for catastrophic energy non-conservation
228  resultNuc = CheckResult(thePro, *nucleus, resultNuc);
229 
230  if(reentryCount>100) {
232  ed << "Call for " << model->GetModelName() << G4endl;
233  ed << "Target element "<<elm->GetName()<<" Z= "
234  << nucleus->GetZ_asInt()
235  << " A= " << nucleus->GetA_asInt() << G4endl;
236  DumpState(track,"ApplyYourself",ed);
237  ed << " ApplyYourself does not completed after 100 attempts" << G4endl;
238  G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006",
239  FatalException, ed);
240  }
241  // Loop checking, 06-Aug-2015, Vladimir Ivanchenko
242  } while(!resultNuc);
243 
244  edep = resultNuc->GetLocalEnergyDeposit();
245  size_t nnuc = resultNuc->GetNumberOfSecondaries();
246 
247  // add delay time of capture
248  for(size_t i=0; i<nnuc; ++i) {
249  G4HadSecondary* sec = resultNuc->GetSecondary(i);
250  sec->SetTime(capTime + sec->GetTime());
251  }
252 
253  nSecondaries += nnuc;
254  result->AddSecondaries(resultNuc);
255  resultNuc->Clear();
256  }
257 
258  // Fill results
259  //
262  theTotalResult->SetNumberOfSecondaries(nSecondaries);
263  G4double w = track.GetWeight();
265  for(G4int i=0; i<nSecondaries; ++i) {
266  G4HadSecondary* sec = result->GetSecondary(i);
267 
268  // add track global time to the reaction time
269  G4double time = sec->GetTime();
270  if(time < 0.0) { time = 0.0; }
271  time += time0;
272 
273  // create secondary track
274  G4Track* t = new G4Track(sec->GetParticle(),
275  time,
276  track.GetPosition());
277  t->SetWeight(w*sec->GetWeight());
278 
279  // use SetCreatorModelIndex to "label" the track
280  if (i<nEmCascadeSec) {
281  t->SetCreatorModelIndex(emcID);
282  } else if (nuclearCapture) {
283  t->SetCreatorModelIndex(ncID);
284  } else {
285  t->SetCreatorModelIndex(dioID);
286  }
287 
290  }
291  result->Clear();
292 
293  if (epReportLevel != 0) {
294  CheckEnergyMomentumConservation(track, *nucleus);
295  }
296  return theTotalResult;
297 }
G4double G4ParticleHPJENDLHEData::G4double result
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
virtual G4Element * SelectZandA(const G4Track &track, G4Nucleus *)
void AddSecondaries(const std::vector< G4HadSecondary > &addSecs)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4HadSecondary * GetSecondary(size_t i)
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
void SetGlobalTime(G4double t)
G4HadProjectile thePro
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
const G4String & GetModelName() const
G4double GetTime() const
int G4int
Definition: G4Types.hh:78
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetCreatorModelIndex(G4int idx)
void ProposeWeight(G4double finalWeight)
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
G4ParticleChange * theTotalResult
bool G4bool
Definition: G4Types.hh:79
G4double GetGlobalTime() const
const G4int n
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
G4Material * GetMaterial() const
void Initialise(const G4Track &aT)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual void Initialize(const G4Track &)
void SetTime(G4double aT)
G4Nucleus * GetTargetNucleusPointer()
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * GetParticle()
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void SetBoundEnergy(G4double e)
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
G4HadronicInteraction * ChooseHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, G4Material *aMaterial, G4Element *anElement)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
const G4String & GetName() const
Definition: G4Element.hh:127
const XML_Char XML_Content * model
Definition: expat.h:151
G4int GetNumberOfSecondaries() const
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4double GetLocalEnergyDeposit() const
G4HadFinalState * CheckResult(const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
G4HadFinalStateStatus GetStatusChange() const
G4double GetWeight() const

Here is the call graph for this function:

G4double G4HadronStoppingProcess::AtRestGetPhysicalInteractionLength ( const G4Track track,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 116 of file G4HadronStoppingProcess.cc.

118 {
119  *condition = NotForced;
120  return 0.0;
121 }
G4double condition(const G4ErrorSymMatrix &m)
void G4HadronStoppingProcess::BuildPhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 109 of file G4HadronStoppingProcess.cc.

110 {
112 }
static G4HadronicProcessStore * Instance()
void PrintInfo(const G4ParticleDefinition *)

Here is the call graph for this function:

G4double G4HadronStoppingProcess::GetMeanLifeTime ( const G4Track ,
G4ForceCondition  
)
inlineprotected

Definition at line 99 of file G4HadronStoppingProcess.hh.

100  { return -1.0; }
G4bool G4HadronStoppingProcess::IsApplicable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4VProcess.

Reimplemented in G4MuonMinusCapture, G4HadronicAbsorptionFritiof, and G4HadronicAbsorptionBertini.

Definition at line 91 of file G4HadronStoppingProcess.cc.

92 {
93  return (p.GetPDGCharge() < 0.0);
94 }
G4double GetPDGCharge() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronStoppingProcess::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 125 of file G4HadronStoppingProcess.cc.

127 {
128  *condition = NotForced;
129  return DBL_MAX;
130 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83
void G4HadronStoppingProcess::PreparePhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 99 of file G4HadronStoppingProcess.cc.

100 {
102  emcID = G4PhysicsModelCatalog::Register(G4String((GetProcessName() + "_EMCascade")));
103  ncID = G4PhysicsModelCatalog::Register(G4String((GetProcessName() + "_NuclearCapture")));
105 }
static G4HadronicProcessStore * Instance()
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
static G4int Register(const G4String &)

Here is the call graph for this function:

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

Reimplemented from G4HadronicProcess.

Reimplemented in G4MuonMinusCapture, G4HadronicAbsorptionFritiof, and G4HadronicAbsorptionBertini.

Definition at line 301 of file G4HadronStoppingProcess.cc.

302 {
303  outFile << "Base process for negatively charged particle capture at rest.\n";
304 }
void G4HadronStoppingProcess::SetBoundDecay ( G4HadronicInteraction ptr)
inline

Definition at line 138 of file G4HadronStoppingProcess.hh.

139 {
140  fBoundDecay = ptr;
141 }

Here is the caller graph for this function:

void G4HadronStoppingProcess::SetElementSelector ( G4ElementSelector ptr)
inline

Definition at line 123 of file G4HadronStoppingProcess.hh.

124 {
125  if(fElementSelector != ptr) {
126  delete fElementSelector;
127  fElementSelector = ptr;
128  }
129 }
void G4HadronStoppingProcess::SetEmCascade ( G4HadronicInteraction ptr)
inline

Definition at line 132 of file G4HadronStoppingProcess.hh.

133 {
134  fEmCascade = ptr;
135 }

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