G4HadronStoppingProcess Class Reference

#include <G4HadronStoppingProcess.hh>

Inheritance diagram for G4HadronStoppingProcess:

G4HadronicProcess G4VDiscreteProcess G4VProcess G4HadronicAbsorptionBertini G4HadronicAbsorptionFritiof G4MuonMinusCapture G4KaonMinusAbsorptionBertini G4PiMinusAbsorptionBertini G4SigmaMinusAbsorptionBertini G4AntiProtonAbsorptionFritiof G4AntiSigmaPlusAbsorptionFritiof

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)

Protected Member Functions

G4double GetMeanLifeTime (const G4Track &, G4ForceCondition *)

Detailed Description

Definition at line 62 of file G4HadronStoppingProcess.hh.


Constructor & Destructor Documentation

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

Definition at line 63 of file G4HadronStoppingProcess.cc.

References G4VProcess::enableAtRestDoIt, G4VProcess::enablePostStepDoIt, G4HadronicProcessStore::Instance(), and G4HadronicProcessStore::RegisterExtraProcess().

00064   : G4HadronicProcess(name, fHadronAtRest)
00065 {
00066   // Modify G4VProcess flags to emulate G4VRest instead of G4VDiscrete
00067   enableAtRestDoIt = true;
00068   enablePostStepDoIt = false;
00069 
00070   fElementSelector = new G4ElementSelector();
00071   fEmCascade = new G4EmCaptureCascade();        // Owned by InteractionRegistry
00072   fBoundDecay = 0;
00073   G4HadronicProcessStore::Instance()->RegisterExtraProcess(this);
00074 }

G4HadronStoppingProcess::~G4HadronStoppingProcess (  )  [virtual]

Definition at line 78 of file G4HadronStoppingProcess.cc.

References G4HadronicProcessStore::DeRegisterExtraProcess(), and G4HadronicProcessStore::Instance().

00079 {
00080   G4HadronicProcessStore::Instance()->DeRegisterExtraProcess(this);
00081   delete fElementSelector;
00082   // NOTE: fEmCascade and fEmBoundDecay owned by registry, not locally
00083 }


Member Function Documentation

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

Reimplemented from G4VDiscreteProcess.

Definition at line 126 of file G4HadronStoppingProcess.cc.

References G4HadFinalState::AddSecondaries(), G4ParticleChange::AddSecondary(), G4HadronicInteraction::ApplyYourself(), G4HadronicProcess::CheckEnergyMomentumConservation(), G4HadronicProcess::CheckResult(), G4HadronicProcess::ChooseHadronicInteraction(), G4HadFinalState::Clear(), G4HadronicProcess::DumpState(), G4HadronicProcess::epReportLevel, FatalException, fStopAndKill, G4endl, G4Exception(), G4Nucleus::GetA_asInt(), G4Track::GetGlobalTime(), G4HadFinalState::GetLocalEnergyDeposit(), G4Track::GetMaterial(), G4HadronicInteraction::GetModelName(), G4Element::GetName(), G4HadFinalState::GetNumberOfSecondaries(), G4Track::GetPosition(), G4HadFinalState::GetSecondary(), G4HadFinalState::GetStatusChange(), G4HadronicProcess::GetTargetNucleusPointer(), G4Track::GetTouchableHandle(), G4Track::GetWeight(), G4Nucleus::GetZ_asInt(), G4HadProjectile::Initialise(), G4ParticleChange::Initialize(), CLHEP::detail::n, G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), G4VParticleChange::ProposeWeight(), G4ElementSelector::SelectZandA(), G4HadProjectile::SetBoundEnergy(), G4VParticleChange::SetNumberOfSecondaries(), G4Track::SetTouchableHandle(), G4Track::SetWeight(), stopAndKill, G4HadronicProcess::thePro, and G4HadronicProcess::theTotalResult.

00128 {
00129   // if primary is not Alive then do nothing
00130   theTotalResult->Initialize(track);
00131 
00132   G4Nucleus* nucleus = GetTargetNucleusPointer();
00133   G4Element* elm = fElementSelector->SelectZandA(track, nucleus);
00134 
00135   G4HadFinalState* result = 0;
00136   thePro.Initialise(track);
00137   G4double time0 = track.GetGlobalTime();
00138   G4bool nuclearCapture = true;
00139 
00140   // Do the electromagnetic cascade in the nuclear field.
00141   // EM cascade should keep G4HadFinalState object,
00142   // because it will not be deleted at the end of this method
00143   //
00144   result = fEmCascade->ApplyYourself(thePro, *nucleus);
00145   G4double ebound = result->GetLocalEnergyDeposit(); 
00146   G4double edep = 0.0; 
00147   G4int nSecondaries = result->GetNumberOfSecondaries();
00148 
00149   // Try decay from bound level 
00150   // For mu- the time of projectile should be changed.
00151   // Decay should keep G4HadFinalState object,
00152   // because it will not be deleted at the end of this method
00153   //
00154   thePro.SetBoundEnergy(ebound);
00155   if(fBoundDecay) {
00156     G4HadFinalState* resultDecay = 
00157       fBoundDecay->ApplyYourself(thePro, *nucleus);
00158     G4int n = resultDecay->GetNumberOfSecondaries();
00159     if(0 < n) {
00160       nSecondaries += n;
00161       result->AddSecondaries(resultDecay);
00162     } 
00163     if(resultDecay->GetStatusChange() == stopAndKill) {
00164       nuclearCapture = false; 
00165     }
00166     resultDecay->Clear();
00167   }
00168 
00169   if(nuclearCapture) {
00170 
00171     // select model
00172     G4HadronicInteraction* model = 0;
00173     try {
00174       model = ChooseHadronicInteraction(0.0, track.GetMaterial(), elm);
00175     }
00176     catch(G4HadronicException & aE) {
00177       G4ExceptionDescription ed;
00178       ed << "Target element "<<elm->GetName()<<"  Z= " 
00179          << nucleus->GetZ_asInt() << "  A= " 
00180          << nucleus->GetA_asInt() << G4endl;
00181       DumpState(track,"ChooseHadronicInteraction",ed);
00182       ed << " No HadronicInteraction found out" << G4endl;
00183       G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had005", 
00184                   FatalException, ed);
00185     }
00186 
00187     G4HadFinalState* resultNuc = 0;
00188     G4int reentryCount = 0;
00189     do {
00190       // sample final state
00191       // nuclear interaction should keep G4HadFinalState object
00192       // model should define time of each secondary particle
00193       try {
00194         resultNuc = model->ApplyYourself(thePro, *nucleus);
00195         ++reentryCount;
00196       }
00197       catch(G4HadronicException aR) {
00198         G4ExceptionDescription ed;
00199         ed << "Call for " << model->GetModelName() << G4endl;
00200         ed << "Target element "<<elm->GetName()<<"  Z= " 
00201            << nucleus->GetZ_asInt() 
00202            << "  A= " << nucleus->GetA_asInt() << G4endl;
00203         DumpState(track,"ApplyYourself",ed);
00204         ed << " ApplyYourself failed" << G4endl;
00205         G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006", 
00206                     FatalException, ed);
00207       }
00208 
00209       // Check the result for catastrophic energy non-conservation
00210       resultNuc = CheckResult(thePro, *nucleus, resultNuc);
00211 
00212       if(reentryCount>100) {
00213         G4ExceptionDescription ed;
00214         ed << "Call for " << model->GetModelName() << G4endl;
00215         ed << "Target element "<<elm->GetName()<<"  Z= " 
00216            << nucleus->GetZ_asInt() 
00217            << "  A= " << nucleus->GetA_asInt() << G4endl;
00218         DumpState(track,"ApplyYourself",ed);
00219         ed << " ApplyYourself does not completed after 100 attempts" << G4endl;
00220         G4Exception("G4HadronStoppingProcess::AtRestDoIt", "had006", 
00221                     FatalException, ed);  
00222       }
00223     }
00224     while(!resultNuc);
00225 
00226     edep = resultNuc->GetLocalEnergyDeposit();
00227     nSecondaries += resultNuc->GetNumberOfSecondaries();
00228     result->AddSecondaries(resultNuc); 
00229     resultNuc->Clear();
00230   }
00231 
00232   // Fill results
00233   //
00234   theTotalResult->ProposeTrackStatus(fStopAndKill);
00235   theTotalResult->ProposeLocalEnergyDeposit(edep);  
00236   theTotalResult->SetNumberOfSecondaries(nSecondaries);
00237   G4double w  = track.GetWeight();
00238   theTotalResult->ProposeWeight(w);
00239   for(G4int i=0; i<nSecondaries; ++i) {
00240     G4HadSecondary* sec = result->GetSecondary(i);
00241 
00242     // add track global time to the reaction time
00243     G4double time = sec->GetTime();
00244     if(time < 0.0) { time = 0.0; }
00245     time += time0;
00246 
00247     // create secondary track
00248     G4Track* t = new G4Track(sec->GetParticle(),
00249                              time, 
00250                              track.GetPosition());
00251     t->SetWeight(w*sec->GetWeight());
00252     t->SetTouchableHandle(track.GetTouchableHandle());
00253     theTotalResult->AddSecondary(t);
00254   }
00255   result->Clear();
00256 
00257   if (epReportLevel != 0) { 
00258     CheckEnergyMomentumConservation(track, *nucleus);
00259   }
00260   return theTotalResult;
00261 }

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

Reimplemented from G4VDiscreteProcess.

Definition at line 108 of file G4HadronStoppingProcess.cc.

References NotForced.

00110 {
00111   *condition = NotForced;
00112   return 0.0;
00113 }

void G4HadronStoppingProcess::BuildPhysicsTable ( const G4ParticleDefinition  )  [virtual]

Reimplemented from G4HadronicProcess.

Definition at line 101 of file G4HadronStoppingProcess.cc.

References G4HadronicProcessStore::Instance(), and G4HadronicProcessStore::PrintInfo().

00102 {
00103   G4HadronicProcessStore::Instance()->PrintInfo(&p);
00104 }

G4double G4HadronStoppingProcess::GetMeanLifeTime ( const G4Track ,
G4ForceCondition  
) [inline, protected]

Definition at line 98 of file G4HadronStoppingProcess.hh.

00099                                                             { return -1.0; }

G4bool G4HadronStoppingProcess::IsApplicable ( const G4ParticleDefinition  )  [virtual]

Reimplemented from G4VProcess.

Reimplemented in G4HadronicAbsorptionBertini, G4HadronicAbsorptionFritiof, and G4MuonMinusCapture.

Definition at line 87 of file G4HadronStoppingProcess.cc.

References G4ParticleDefinition::GetPDGCharge().

Referenced by G4HadronicAbsorptionBertini::IsApplicable().

00088 {
00089   return (p.GetPDGCharge() < 0.0);
00090 }

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

Reimplemented from G4VDiscreteProcess.

Definition at line 117 of file G4HadronStoppingProcess.cc.

References DBL_MAX, and NotForced.

00119 {
00120   *condition = NotForced;
00121   return DBL_MAX;
00122 }

void G4HadronStoppingProcess::PreparePhysicsTable ( const G4ParticleDefinition  )  [virtual]

Reimplemented from G4HadronicProcess.

Definition at line 94 of file G4HadronStoppingProcess.cc.

References G4HadronicProcessStore::Instance(), and G4HadronicProcessStore::RegisterParticleForExtraProcess().

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

Reimplemented from G4HadronicProcess.

Reimplemented in G4HadronicAbsorptionBertini, G4HadronicAbsorptionFritiof, and G4MuonMinusCapture.

Definition at line 265 of file G4HadronStoppingProcess.cc.

00266 {
00267   outFile << "Base process for negatively charged particle capture at rest.\n";
00268 }

void G4HadronStoppingProcess::SetBoundDecay ( G4HadronicInteraction ptr  )  [inline]

Definition at line 133 of file G4HadronStoppingProcess.hh.

Referenced by G4MuonMinusCapture::G4MuonMinusCapture().

00134 {
00135   fBoundDecay = ptr;
00136 }

void G4HadronStoppingProcess::SetElementSelector ( G4ElementSelector ptr  )  [inline]

Definition at line 118 of file G4HadronStoppingProcess.hh.

00119 {
00120   if(fElementSelector != ptr) {
00121     delete fElementSelector;
00122     fElementSelector = ptr;
00123   }
00124 }

void G4HadronStoppingProcess::SetEmCascade ( G4HadronicInteraction ptr  )  [inline]

Definition at line 127 of file G4HadronStoppingProcess.hh.

00128 {
00129   fEmCascade = ptr;
00130 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:09 2013 for Geant4 by  doxygen 1.4.7