#include <G4NuclearStopping.hh>
Inheritance diagram for G4NuclearStopping:
Public Member Functions | |
G4NuclearStopping (const G4String &processName="nuclearStopping") | |
virtual | ~G4NuclearStopping () |
G4bool | IsApplicable (const G4ParticleDefinition &p) |
G4double | AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection) |
G4VParticleChange * | AlongStepDoIt (const G4Track &track, const G4Step &step) |
void | PrintInfo () |
Protected Member Functions | |
void | InitialiseProcess (const G4ParticleDefinition *) |
Definition at line 65 of file G4NuclearStopping.hh.
G4NuclearStopping::G4NuclearStopping | ( | const G4String & | processName = "nuclearStopping" |
) |
Definition at line 54 of file G4NuclearStopping.cc.
References G4VProcess::enableAlongStepDoIt, G4VProcess::enablePostStepDoIt, fNuclearStopping, G4VEmProcess::SetBuildTableFlag(), and G4VProcess::SetProcessSubType().
00055 : G4VEmProcess(processName) 00056 { 00057 isInitialized = false; 00058 SetProcessSubType(fNuclearStopping); 00059 SetBuildTableFlag(false); 00060 enableAlongStepDoIt = true; 00061 enablePostStepDoIt = false; 00062 modelICRU49 = 0; 00063 }
G4NuclearStopping::~G4NuclearStopping | ( | ) | [virtual] |
G4VParticleChange * G4NuclearStopping::AlongStepDoIt | ( | const G4Track & | track, | |
const G4Step & | step | |||
) | [virtual] |
Reimplemented from G4VDiscreteProcess.
Definition at line 105 of file G4NuclearStopping.cc.
References G4VEmModel::ComputeDEDXPerVolume(), G4MaterialCutsCouple::GetIndex(), G4StepPoint::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4Track::GetMaterialCutsCouple(), G4Track::GetParticleDefinition(), G4ParticleDefinition::GetPDGCharge(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4Step::GetStepLength(), G4ParticleChangeForLoss::InitializeForAlongStep(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeNonIonizingEnergyDeposit(), G4VEmProcess::SelectModel(), G4ICRU49NuclearStoppingModel::SetFluctuationFlag(), and G4ParticleChangeForLoss::SetProposedKineticEnergy().
00107 { 00108 nParticleChange.InitializeForAlongStep(track); 00109 G4double T2 = step.GetPostStepPoint()->GetKineticEnergy(); 00110 00111 const G4ParticleDefinition* part = track.GetParticleDefinition(); 00112 G4double Z = std::fabs(part->GetPDGCharge()/eplus); 00113 00114 if(T2 > 0.0 && T2*proton_mass_c2 < Z*Z*MeV*part->GetPDGMass()) { 00115 00116 G4double length = step.GetStepLength(); 00117 if(length > 0.0) { 00118 00119 // primary 00120 G4double T1= step.GetPreStepPoint()->GetKineticEnergy(); 00121 G4double T = 0.5*(T1 + T2); 00122 const G4MaterialCutsCouple* couple = track.GetMaterialCutsCouple(); 00123 G4VEmModel* mod = SelectModel(T, couple->GetIndex()); 00124 00125 // sample stopping 00126 if(modelICRU49) { modelICRU49->SetFluctuationFlag(true); } 00127 G4double nloss = 00128 length*mod->ComputeDEDXPerVolume(couple->GetMaterial(), part, T); 00129 if(nloss > T1) { nloss = T1; } 00130 nParticleChange.SetProposedKineticEnergy(T1 - nloss); 00131 nParticleChange.ProposeLocalEnergyDeposit(nloss); 00132 nParticleChange.ProposeNonIonizingEnergyDeposit(nloss); 00133 if(modelICRU49) { modelICRU49->SetFluctuationFlag(false); } 00134 } 00135 } 00136 return &nParticleChange; 00137 }
G4double G4NuclearStopping::AlongStepGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4double | currentMinimumStep, | |||
G4double & | proposedSafety, | |||
G4GPILSelection * | selection | |||
) | [virtual] |
Reimplemented from G4VDiscreteProcess.
Definition at line 96 of file G4NuclearStopping.cc.
References DBL_MAX, and NotCandidateForSelection.
00098 { 00099 *selection = NotCandidateForSelection; 00100 return DBL_MAX; 00101 }
void G4NuclearStopping::InitialiseProcess | ( | const G4ParticleDefinition * | ) | [protected, virtual] |
Implements G4VEmProcess.
Definition at line 79 of file G4NuclearStopping.cc.
References G4VEmProcess::AddEmModel(), G4VEmProcess::EmModel(), G4VEmProcess::SetEmModel(), and G4VEmModel::SetParticleChange().
00080 { 00081 if(!isInitialized) { 00082 isInitialized = true; 00083 00084 if(!EmModel(1)) { 00085 modelICRU49 = new G4ICRU49NuclearStoppingModel(); 00086 SetEmModel(modelICRU49); 00087 } 00088 AddEmModel(1, EmModel()); 00089 EmModel()->SetParticleChange(&nParticleChange); 00090 } 00091 }
G4bool G4NuclearStopping::IsApplicable | ( | const G4ParticleDefinition & | p | ) | [virtual] |
Implements G4VEmProcess.
Definition at line 72 of file G4NuclearStopping.cc.
References G4ParticleDefinition::GetPDGCharge(), and G4ParticleDefinition::IsShortLived().
00073 { 00074 return (p.GetPDGCharge() != 0.0 && !p.IsShortLived()); 00075 }
void G4NuclearStopping::PrintInfo | ( | ) | [virtual] |