Geant4-11
Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes
G4NeutronElectronElModel Class Reference

#include <G4NeutronElectronElModel.hh>

Inheritance diagram for G4NeutronElectronElModel:
G4HadronElastic G4HadronicInteraction

Public Member Functions

void ActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Material *aMaterial)
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4double CalculateAm (G4double momentum)
 
G4double ComputeMomentumCMS (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
void DeActivateFor (const G4Element *anElement)
 
void DeActivateFor (const G4Material *aMaterial)
 
 G4NeutronElectronElModel (const G4String &name="n-e-elastic")
 
G4double GetAm ()
 
G4double GetCutEnergy ()
 
virtual std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
virtual const std::pair< G4double, G4doubleGetFatalEnergyCheckLevels () const
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
const G4StringGetModelName () const
 
G4double GetRecoilEnergyThreshold () const
 
G4double GetSlopeCof (const G4int pdg)
 
G4double GetTransfer (G4int iTkin, G4int iTransfer, G4double position)
 
G4int GetVerboseLevel () const
 
void Initialise ()
 
virtual void InitialiseModel ()
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4bool IsBlocked (const G4Element *anElement) const
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4double LowestEnergyLimit () const
 
virtual void ModelDescription (std::ostream &) const
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A) override
 
G4double SampleSin2HalfTheta (G4double Tkin)
 
void SetCutEnergy (G4double ec)
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
void SetLowestEnergyLimit (G4double value)
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
void SetRecoilEnergyThreshold (G4double val)
 
void SetVerboseLevel (G4int value)
 
G4double XscIntegrand (G4double x)
 
virtual ~G4NeutronElectronElModel ()
 

Protected Member Functions

void Block ()
 
G4bool IsBlocked () const
 
void SetModelName (const G4String &nam)
 

Protected Attributes

G4bool isBlocked
 
G4double pLocalTmax
 
G4int secID
 
G4double theMaxEnergy
 
G4double theMinEnergy
 
G4HadFinalState theParticleChange
 
G4int verboseLevel
 

Private Attributes

std::pair< G4double, G4doubleepCheckLevels
 
G4double fAm
 
G4int fAngleBin
 
G4PhysicsTablefAngleTable
 
G4double fCutEnergy
 
G4double fee
 
G4double fee2
 
G4int fEnergyBin
 
G4PhysicsLogVectorfEnergyVector
 
G4double fM
 
G4double fM2
 
G4double fMaxEnergy
 
G4double fme
 
G4double fme2
 
G4double fMinEnergy
 
G4double fMv2
 
G4double lowestEnergyLimit
 
G4int nwarn
 
G4double recoilEnergyThreshold
 
G4HadronicInteractionRegistryregistry
 
G4ParticleDefinitiontheAlpha
 
std::vector< const G4Material * > theBlockedList
 
std::vector< const G4Element * > theBlockedListElements
 
G4ParticleDefinitiontheDeuteron
 
G4ParticleDefinitiontheElectron
 
std::vector< std::pair< G4double, const G4Material * > > theMaxEnergyList
 
std::vector< std::pair< G4double, const G4Element * > > theMaxEnergyListElements
 
std::vector< std::pair< G4double, const G4Material * > > theMinEnergyList
 
std::vector< std::pair< G4double, const G4Element * > > theMinEnergyListElements
 
G4String theModelName
 
G4ParticleDefinitiontheNeutron
 
G4ParticleDefinitiontheProton
 

Detailed Description

Definition at line 51 of file G4NeutronElectronElModel.hh.

Constructor & Destructor Documentation

◆ G4NeutronElectronElModel()

G4NeutronElectronElModel::G4NeutronElectronElModel ( const G4String name = "n-e-elastic")

Definition at line 49 of file G4NeutronElectronElModel.cc.

51{
53
54 // neutron magneton squared
55
56 fM = neutron_mass_c2; // neutron mass
57 fM2 = fM*fM;
59 fme2 = fme*fme;
60 fMv2 = 0.7056*GeV*GeV;
61
62 SetMinEnergy( 0.001*GeV );
63 SetMaxEnergy( 10.*TeV );
65
67 // PDG2016: sin^2 theta Weinberg
68
69 fEnergyBin = 200;
70 fMinEnergy = 1.*MeV;
71 fMaxEnergy = 10000.*GeV;
73
74 fAngleBin = 500;
75 fAngleTable = 0;
76
77 fCutEnergy = 0.; // default value
78
79 Initialise();
80}
static constexpr double eV
Definition: G4SIunits.hh:201
static constexpr double GeV
Definition: G4SIunits.hh:203
static constexpr double MeV
Definition: G4SIunits.hh:200
static constexpr double TeV
Definition: G4SIunits.hh:204
static G4Electron * Electron()
Definition: G4Electron.cc:93
G4HadronElastic(const G4String &name="hElasticLHEP")
void SetLowestEnergyLimit(G4double value)
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
G4PhysicsLogVector * fEnergyVector
G4ParticleDefinition * theElectron
static G4int GetModelID(const G4int modelIndex)
const char * name(G4int ptype)
float electron_mass_c2
Definition: hepunit.py:273
float neutron_mass_c2
Definition: hepunit.py:275

References G4Electron::Electron(), source.hepunit::electron_mass_c2, eV, fAngleBin, fAngleTable, fCutEnergy, fEnergyBin, fEnergyVector, fM, fM2, fMaxEnergy, fme, fme2, fMinEnergy, fMv2, G4PhysicsModelCatalog::GetModelID(), GeV, Initialise(), MeV, G4InuclParticleNames::name(), source.hepunit::neutron_mass_c2, G4HadronElastic::secID, G4HadronElastic::SetLowestEnergyLimit(), G4HadronicInteraction::SetMaxEnergy(), G4HadronicInteraction::SetMinEnergy(), TeV, and theElectron.

Referenced by Initialise().

◆ ~G4NeutronElectronElModel()

G4NeutronElectronElModel::~G4NeutronElectronElModel ( )
virtual

Definition at line 84 of file G4NeutronElectronElModel.cc.

85{
86 if( fEnergyVector )
87 {
88 delete fEnergyVector;
89 fEnergyVector = nullptr;
90 }
91 if( fAngleTable )
92 {
94 delete fAngleTable;
95 fAngleTable = nullptr;
96 }
97}
void clearAndDestroy()

References G4PhysicsTable::clearAndDestroy(), fAngleTable, and fEnergyVector.

Member Function Documentation

◆ ActivateFor() [1/2]

void G4HadronicInteraction::ActivateFor ( const G4Element anElement)
inlineinherited

◆ ActivateFor() [2/2]

void G4HadronicInteraction::ActivateFor ( const G4Material aMaterial)
inlineinherited

◆ ApplyYourself()

G4HadFinalState * G4NeutronElectronElModel::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Reimplemented from G4HadronElastic.

Definition at line 262 of file G4NeutronElectronElModel.cc.

264{
266
267 const G4HadProjectile* aParticle = &aTrack;
268 G4double Tkin = aParticle->GetKineticEnergy();
269 fAm = CalculateAm( Tkin);
270 // G4double En = aParticle->GetTotalEnergy();
271
272 if( Tkin <= LowestEnergyLimit() )
273 {
276 return &theParticleChange;
277 }
278 // sample e-scattering angle and make final state in lab frame
279
280 G4double sin2ht = SampleSin2HalfTheta( Tkin); // in n-rrest frame
281
282 // G4cout<<"sin2ht = "<<sin2ht<<G4endl;
283
284 G4double eTkin = fee; // fM;
285
286 eTkin /= 1.+2.*fee*sin2ht/fM; // fme/En + 2*sin2ht;
287
288 eTkin -= fme;
289
290 // G4cout<<"eTkin = "<<eTkin<<G4endl;
291
292 if( eTkin > fCutEnergy )
293 {
294 G4double ePlab = sqrt( eTkin*(eTkin + 2.*fme) );
295
296 // G4cout<<"ePlab = "<<ePlab<<G4endl;
297
298 G4double cost = 1. - 2*sin2ht;
299
300 if( cost > 1. ) cost = 1.;
301 if( cost < -1. ) cost = -1.;
302
303 G4double sint = std::sqrt( (1.0 - cost)*(1.0 + cost) );
305
306 G4ThreeVector eP( sint*std::cos(phi), sint*std::sin(phi), cost );
307 eP *= ePlab;
308 G4LorentzVector lvt2( eP, eTkin + electron_mass_c2 ); // recoil e- in n-rest frame
309
310 G4LorentzVector lvp1 = aParticle->Get4Momentum();
311 G4LorentzVector lvt1(0.,0.,0.,electron_mass_c2);
312 G4LorentzVector lvsum = lvp1+lvt1;
313
314 G4ThreeVector bst = lvp1.boostVector();
315 lvt2.boost(bst);
316
317 // G4cout<<"lvt2 = "<<lvt2<<G4endl;
318
319 G4DynamicParticle * aSec = new G4DynamicParticle( theElectron, lvt2 );
321
322 G4LorentzVector lvp2 = lvsum-lvt2;
323
324 // G4cout<<"lvp2 = "<<lvp2<<G4endl;
325
326 G4double Tkin2 = lvp2.e()-aParticle->GetDefinition()->GetPDGMass();
329 }
330 else if( eTkin > 0.0 )
331 {
333 Tkin -= eTkin;
334
335 if( Tkin > 0. )
336 {
339 }
340 }
341 else
342 {
345 }
346 return &theParticleChange;
347}
double G4double
Definition: G4Types.hh:83
#define G4UniformRand()
Definition: Randomize.hh:52
Hep3Vector unit() const
Hep3Vector boostVector() const
Hep3Vector vect() const
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
void SetEnergyChange(G4double anEnergy)
void SetMomentumChange(const G4ThreeVector &aV)
void SetLocalEnergyDeposit(G4double aE)
const G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
G4double LowestEnergyLimit() const
G4double SampleSin2HalfTheta(G4double Tkin)
G4double CalculateAm(G4double momentum)
static constexpr double twopi
Definition: SystemOfUnits.h:56

References G4HadFinalState::AddSecondary(), CLHEP::HepLorentzVector::boost(), CLHEP::HepLorentzVector::boostVector(), CalculateAm(), G4HadFinalState::Clear(), CLHEP::HepLorentzVector::e(), source.hepunit::electron_mass_c2, fAm, fCutEnergy, fee, fM, fme, G4UniformRand, G4HadProjectile::Get4Momentum(), G4HadProjectile::GetDefinition(), G4HadProjectile::GetKineticEnergy(), G4ParticleDefinition::GetPDGMass(), G4HadronElastic::LowestEnergyLimit(), SampleSin2HalfTheta(), G4HadronElastic::secID, G4HadFinalState::SetEnergyChange(), G4HadFinalState::SetLocalEnergyDeposit(), G4HadFinalState::SetMomentumChange(), theElectron, G4HadronicInteraction::theParticleChange, CLHEP::twopi, CLHEP::Hep3Vector::unit(), and CLHEP::HepLorentzVector::vect().

◆ Block()

void G4HadronicInteraction::Block ( )
inlineprotectedinherited

◆ BuildPhysicsTable()

void G4HadronicInteraction::BuildPhysicsTable ( const G4ParticleDefinition )
virtualinherited

◆ CalculateAm()

G4double G4NeutronElectronElModel::CalculateAm ( G4double  momentum)
inline

Definition at line 101 of file G4NeutronElectronElModel.hh.

102{
103 fee = (Tkin+fM)*fme/fM;
104 // G4cout<<"fee = "<<fee<<" MeV"<<G4endl;
105 fee2 = fee*fee;
106 G4double momentum = std::sqrt( fee2 - fme2 );
107 G4double k = momentum/CLHEP::hbarc;
108 G4double ch = 1.13;
109 G4double zn = 1.77*k*CLHEP::Bohr_radius;
110 G4double zn2 = zn*zn;
111 fAm = ch/zn2;
112
113 return fAm;
114}
static constexpr double Bohr_radius
static constexpr double hbarc

References CLHEP::Bohr_radius, fAm, fee, fee2, fM, fme, fme2, and CLHEP::hbarc.

Referenced by ApplyYourself(), and Initialise().

◆ ComputeMomentumCMS()

G4double G4HadronElastic::ComputeMomentumCMS ( const G4ParticleDefinition p,
G4double  plab,
G4int  Z,
G4int  A 
)
inlineinherited

Definition at line 102 of file G4HadronElastic.hh.

104{
105 G4double m1 = p->GetPDGMass();
106 G4double m12= m1*m1;
108 return plab*mass2/std::sqrt(m12 + mass2*mass2 + 2.*mass2*std::sqrt(m12 + plab*plab));
109}
const G4int Z[17]
const G4double A[17]
static G4double GetNuclearMass(const G4double A, const G4double Z)

References A, G4NucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetPDGMass(), and Z.

◆ DeActivateFor() [1/2]

void G4HadronicInteraction::DeActivateFor ( const G4Element anElement)
inherited

Definition at line 186 of file G4HadronicInteraction.cc.

187{
188 Block();
189 theBlockedListElements.push_back(anElement);
190}
std::vector< const G4Element * > theBlockedListElements

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theBlockedListElements.

◆ DeActivateFor() [2/2]

void G4HadronicInteraction::DeActivateFor ( const G4Material aMaterial)
inherited

Definition at line 180 of file G4HadronicInteraction.cc.

181{
182 Block();
183 theBlockedList.push_back(aMaterial);
184}
std::vector< const G4Material * > theBlockedList

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theBlockedList.

Referenced by G4HadronHElasticPhysics::ConstructProcess().

◆ GetAm()

G4double G4NeutronElectronElModel::GetAm ( )
inline

Definition at line 71 of file G4NeutronElectronElModel.hh.

71{return fAm;};

References fAm.

◆ GetCutEnergy()

G4double G4NeutronElectronElModel::GetCutEnergy ( )
inline

Definition at line 80 of file G4NeutronElectronElModel.hh.

80{return fCutEnergy;};

References fCutEnergy.

◆ GetEnergyMomentumCheckLevels()

std::pair< G4double, G4double > G4HadronicInteraction::GetEnergyMomentumCheckLevels ( ) const
virtualinherited

◆ GetFatalEnergyCheckLevels()

const std::pair< G4double, G4double > G4HadronicInteraction::GetFatalEnergyCheckLevels ( ) const
virtualinherited

Reimplemented in G4FissLib, G4LFission, G4LENDFission, G4ParticleHPCapture, G4ParticleHPElastic, G4ParticleHPFission, G4ParticleHPInelastic, and G4ParticleHPThermalScattering.

Definition at line 210 of file G4HadronicInteraction.cc.

211{
212 // default level of Check
213 return std::pair<G4double, G4double>(2.*perCent, 1. * GeV);
214}
static constexpr double perCent
Definition: G4SIunits.hh:325

References GeV, and perCent.

Referenced by G4HadronicProcess::CheckResult().

◆ GetMaxEnergy() [1/2]

G4double G4HadronicInteraction::GetMaxEnergy ( ) const
inlineinherited

◆ GetMaxEnergy() [2/2]

G4double G4HadronicInteraction::GetMaxEnergy ( const G4Material aMaterial,
const G4Element anElement 
) const
inherited

Definition at line 131 of file G4HadronicInteraction.cc.

133{
134 if(!IsBlocked()) { return theMaxEnergy; }
135 if( IsBlocked(aMaterial) || IsBlocked(anElement) ) { return 0.0; }
136 if(!theMaxEnergyListElements.empty()) {
137 for(auto const& elmlist : theMaxEnergyListElements) {
138 if( anElement == elmlist.second )
139 { return elmlist.first; }
140 }
141 }
142 if(!theMaxEnergyList.empty()) {
143 for(auto const& matlist : theMaxEnergyList) {
144 if( aMaterial == matlist.second )
145 { return matlist.first; }
146 }
147 }
148 return theMaxEnergy;
149}
std::vector< std::pair< G4double, const G4Material * > > theMaxEnergyList
std::vector< std::pair< G4double, const G4Element * > > theMaxEnergyListElements

References G4HadronicInteraction::IsBlocked(), G4HadronicInteraction::theMaxEnergy, G4HadronicInteraction::theMaxEnergyList, and G4HadronicInteraction::theMaxEnergyListElements.

◆ GetMinEnergy() [1/2]

G4double G4HadronicInteraction::GetMinEnergy ( ) const
inlineinherited

◆ GetMinEnergy() [2/2]

G4double G4HadronicInteraction::GetMinEnergy ( const G4Material aMaterial,
const G4Element anElement 
) const
inherited

Definition at line 81 of file G4HadronicInteraction.cc.

83{
84 if(!IsBlocked()) { return theMinEnergy; }
85 if( IsBlocked(aMaterial) || IsBlocked(anElement) ) { return DBL_MAX; }
86 if(!theMinEnergyListElements.empty()) {
87 for(auto const& elmlist : theMinEnergyListElements) {
88 if( anElement == elmlist.second )
89 { return elmlist.first; }
90 }
91 }
92 if(!theMinEnergyList.empty()) {
93 for(auto const & matlist : theMinEnergyList) {
94 if( aMaterial == matlist.second )
95 { return matlist.first; }
96 }
97 }
98 return theMinEnergy;
99}
std::vector< std::pair< G4double, const G4Element * > > theMinEnergyListElements
std::vector< std::pair< G4double, const G4Material * > > theMinEnergyList
#define DBL_MAX
Definition: templates.hh:62

References DBL_MAX, G4HadronicInteraction::IsBlocked(), G4HadronicInteraction::theMinEnergy, G4HadronicInteraction::theMinEnergyList, and G4HadronicInteraction::theMinEnergyListElements.

◆ GetModelName()

const G4String & G4HadronicInteraction::GetModelName ( ) const
inlineinherited

Definition at line 115 of file G4HadronicInteraction.hh.

116 { return theModelName; }

References G4HadronicInteraction::theModelName.

Referenced by G4MuMinusCapturePrecompound::ApplyYourself(), G4HadronElastic::ApplyYourself(), G4INCLXXInterface::ApplyYourself(), G4TheoFSGenerator::ApplyYourself(), G4HadronStoppingProcess::AtRestDoIt(), G4VHadronPhysics::BuildModel(), G4HadronicProcess::CheckEnergyMomentumConservation(), G4HadronicProcess::CheckResult(), G4ChargeExchangePhysics::ConstructProcess(), G4MuonicAtomDecay::DecayIt(), G4LENDModel::DumpLENDTargetInfo(), G4AblaInterface::G4AblaInterface(), G4ElectroVDNuclearModel::G4ElectroVDNuclearModel(), G4EMDissociation::G4EMDissociation(), G4ExcitedStringDecay::G4ExcitedStringDecay(), G4LEHadronProtonElastic::G4LEHadronProtonElastic(), G4LENDModel::G4LENDModel(), G4LENDorBERTModel::G4LENDorBERTModel(), G4LEnp::G4LEnp(), G4LEpp::G4LEpp(), G4LFission::G4LFission(), G4LowEGammaNuclearModel::G4LowEGammaNuclearModel(), G4LowEIonFragmentation::G4LowEIonFragmentation(), G4MuonVDNuclearModel::G4MuonVDNuclearModel(), G4NeutrinoElectronCcModel::G4NeutrinoElectronCcModel(), G4NeutrinoNucleusModel::G4NeutrinoNucleusModel(), G4WilsonAbrasionModel::G4WilsonAbrasionModel(), G4INCLXXInterface::GetDeExcitationModelName(), G4EnergyRangeManager::GetHadronicInteraction(), G4VHighEnergyGenerator::GetProjectileNucleus(), G4NeutronRadCapture::InitialiseModel(), G4BinaryCascade::ModelDescription(), G4LMsdGenerator::ModelDescription(), G4VPartonStringModel::ModelDescription(), G4TheoFSGenerator::ModelDescription(), G4VHadronPhysics::NewModel(), G4NeutrinoElectronProcess::PostStepDoIt(), G4HadronicProcess::PostStepDoIt(), G4ElNeutrinoNucleusProcess::PostStepDoIt(), G4HadronElasticProcess::PostStepDoIt(), G4MuNeutrinoNucleusProcess::PostStepDoIt(), G4HadronicProcessStore::PrintModelHtml(), G4BinaryCascade::PropagateModelDescription(), G4HadronicProcessStore::RegisterInteraction(), and G4LENDModel::returnUnchanged().

◆ GetRecoilEnergyThreshold()

G4double G4HadronicInteraction::GetRecoilEnergyThreshold ( ) const
inlineinherited

◆ GetSlopeCof()

G4double G4HadronElastic::GetSlopeCof ( const G4int  pdg)
inherited

Definition at line 277 of file G4HadronElastic.cc.

278{
279 // The input parameter "pdg" should be the absolute value of the PDG code
280 // (i.e. the same value for a particle and its antiparticle).
281
282 G4double coeff = 1.0;
283
284 // heavy barions
285
286 static const G4double lBarCof1S = 0.88;
287 static const G4double lBarCof2S = 0.76;
288 static const G4double lBarCof3S = 0.64;
289 static const G4double lBarCof1C = 0.784378;
290 static const G4double lBarCofSC = 0.664378;
291 static const G4double lBarCof2SC = 0.544378;
292 static const G4double lBarCof1B = 0.740659;
293 static const G4double lBarCofSB = 0.620659;
294 static const G4double lBarCof2SB = 0.500659;
295
296 if( pdg == 3122 || pdg == 3222 || pdg == 3112 || pdg == 3212 )
297 {
298 coeff = lBarCof1S; // Lambda, Sigma+, Sigma-, Sigma0
299
300 } else if( pdg == 3322 || pdg == 3312 )
301 {
302 coeff = lBarCof2S; // Xi-, Xi0
303 }
304 else if( pdg == 3324)
305 {
306 coeff = lBarCof3S; // Omega
307 }
308 else if( pdg == 4122 || pdg == 4212 || pdg == 4222 || pdg == 4112 )
309 {
310 coeff = lBarCof1C; // LambdaC+, SigmaC+, SigmaC++, SigmaC0
311 }
312 else if( pdg == 4332 )
313 {
314 coeff = lBarCof2SC; // OmegaC
315 }
316 else if( pdg == 4232 || pdg == 4132 )
317 {
318 coeff = lBarCofSC; // XiC+, XiC0
319 }
320 else if( pdg == 5122 || pdg == 5222 || pdg == 5112 || pdg == 5212 )
321 {
322 coeff = lBarCof1B; // LambdaB, SigmaB+, SigmaB-, SigmaB0
323 }
324 else if( pdg == 5332 )
325 {
326 coeff = lBarCof2SB; // OmegaB-
327 }
328 else if( pdg == 5132 || pdg == 5232 ) // XiB-, XiB0
329 {
330 coeff = lBarCofSB;
331 }
332 // heavy mesons Kaons?
333 static const G4double lMesCof1S = 0.82; // Kp/piP kaons?
334 static const G4double llMesCof1C = 0.676568;
335 static const G4double llMesCof1B = 0.610989;
336 static const G4double llMesCof2C = 0.353135;
337 static const G4double llMesCof2B = 0.221978;
338 static const G4double llMesCofSC = 0.496568;
339 static const G4double llMesCofSB = 0.430989;
340 static const G4double llMesCofCB = 0.287557;
341 static const G4double llMesCofEtaP = 0.88;
342 static const G4double llMesCofEta = 0.76;
343
344 if( pdg == 321 || pdg == 311 || pdg == 310 )
345 {
346 coeff = lMesCof1S; //K+-0
347 }
348 else if( pdg == 511 || pdg == 521 )
349 {
350 coeff = llMesCof1B; // BMeson0, BMeson+
351 }
352 else if(pdg == 421 || pdg == 411 )
353 {
354 coeff = llMesCof1C; // DMeson+, DMeson0
355 }
356 else if( pdg == 531 )
357 {
358 coeff = llMesCofSB; // BSMeson0
359 }
360 else if( pdg == 541 )
361 {
362 coeff = llMesCofCB; // BCMeson+-
363 }
364 else if(pdg == 431 )
365 {
366 coeff = llMesCofSC; // DSMeson+-
367 }
368 else if(pdg == 441 || pdg == 443 )
369 {
370 coeff = llMesCof2C; // Etac, JPsi
371 }
372 else if(pdg == 553 )
373 {
374 coeff = llMesCof2B; // Upsilon
375 }
376 else if(pdg == 221 )
377 {
378 coeff = llMesCofEta; // Eta
379 }
380 else if(pdg == 331 )
381 {
382 coeff = llMesCofEtaP; // Eta'
383 }
384 return coeff;
385}

◆ GetTransfer()

G4double G4NeutronElectronElModel::GetTransfer ( G4int  iTkin,
G4int  iTransfer,
G4double  position 
)

Definition at line 194 of file G4NeutronElectronElModel.cc.

195{
196 G4double x1, x2, y1, y2, randTransfer, delta, mean, epsilon = 1.e-6;
197
198 if( iTransfer == 0 || iTransfer == fAngleBin-1 )
199 {
200 randTransfer = (*fAngleTable)(iTkin)->Energy(iTransfer);
201 }
202 else
203 {
204 if ( iTransfer >= G4int((*fAngleTable)(iTkin)->GetVectorLength()) )
205 {
206 iTransfer = (*fAngleTable)(iTkin)->GetVectorLength() - 1;
207 }
208 y1 = (*(*fAngleTable)(iTkin))(iTransfer-1);
209 y2 = (*(*fAngleTable)(iTkin))(iTransfer);
210
211 x1 = (*fAngleTable)(iTkin)->Energy(iTransfer-1);
212 x2 = (*fAngleTable)(iTkin)->Energy(iTransfer);
213
214 delta = y2 - y1;
215 mean = y2 + y1;
216
217 if ( x1 == x2 ) randTransfer = x2;
218 else
219 {
220 if ( delta < epsilon*mean )
221 {
222 randTransfer = x1 + ( x2 - x1 )*G4UniformRand();
223 }
224 else
225 {
226 randTransfer = x1 + ( position - y1 )*( x2 - x1 )/delta; // ( y2 - y1 );
227 }
228 }
229 }
230 return randTransfer;
231}
G4double epsilon(G4double density, G4double temperature)
int G4int
Definition: G4Types.hh:85

References epsilon(), fAngleBin, fAngleTable, and G4UniformRand.

Referenced by SampleSin2HalfTheta().

◆ GetVerboseLevel()

G4int G4HadronicInteraction::GetVerboseLevel ( ) const
inlineinherited

Definition at line 109 of file G4HadronicInteraction.hh.

References G4HadronicInteraction::verboseLevel.

◆ Initialise()

void G4NeutronElectronElModel::Initialise ( )

Definition at line 120 of file G4NeutronElectronElModel.cc.

121{
122 G4double result = 0., sum, Tkin, dt, t1, t2;
123 G4int iTkin, jTransfer;
125
127
128 for( iTkin = 0; iTkin < fEnergyBin; iTkin++)
129 {
130 Tkin = fEnergyVector->GetLowEdgeEnergy(iTkin);
131 fAm = CalculateAm(Tkin);
132 dt = 1./fAngleBin;
133
135
136 sum = 0.;
137
138 for( jTransfer = 0; jTransfer < fAngleBin; jTransfer++)
139 {
140 t1 = dt*jTransfer;
141 t2 = t1 + dt;
142
143 result = integral.Legendre96( this, &G4NeutronElectronElModel::XscIntegrand, t1, t2 );
144
145 sum += result;
146 // G4cout<<sum<<", ";
147 vectorT->PutValue(jTransfer, t1, sum);
148 }
149 // G4cout<<G4endl;
150 fAngleTable->insertAt(iTkin,vectorT);
151 }
152 return;
153}
G4NeutronElectronElModel(const G4String &name="n-e-elastic")
void PutValue(const std::size_t index, const G4double e, const G4double value)
void insertAt(std::size_t, G4PhysicsVector *)
G4double GetLowEdgeEnergy(const std::size_t index) const

References CalculateAm(), fAm, fAngleBin, fAngleTable, fEnergyBin, fEnergyVector, G4NeutronElectronElModel(), G4PhysicsVector::GetLowEdgeEnergy(), G4PhysicsTable::insertAt(), G4PhysicsFreeVector::PutValue(), and XscIntegrand().

Referenced by G4NeutronElectronElModel().

◆ InitialiseModel()

void G4HadronicInteraction::InitialiseModel ( )
virtualinherited

◆ IsApplicable()

G4bool G4NeutronElectronElModel::IsApplicable ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 110 of file G4NeutronElectronElModel.cc.

111{
112 G4String pName = aTrack.GetDefinition()->GetParticleName();
113 G4double energy = aTrack.GetTotalEnergy();
114
115 return (pName == "neutron" && energy >= fMinEnergy && energy <= fMaxEnergy);
116}
G4double GetTotalEnergy() const
const G4String & GetParticleName() const
G4double energy(const ThreeVector &p, const G4double m)

References G4INCL::KinematicsUtils::energy(), fMaxEnergy, fMinEnergy, G4HadProjectile::GetDefinition(), G4ParticleDefinition::GetParticleName(), and G4HadProjectile::GetTotalEnergy().

◆ IsBlocked() [1/3]

G4bool G4HadronicInteraction::IsBlocked ( ) const
inlineprotectedinherited

◆ IsBlocked() [2/3]

G4bool G4HadronicInteraction::IsBlocked ( const G4Element anElement) const
inherited

Definition at line 202 of file G4HadronicInteraction.cc.

203{
204 for (auto const& elm : theBlockedListElements) {
205 if (anElement == elm) return true;
206 }
207 return false;
208}

References G4HadronicInteraction::theBlockedListElements.

◆ IsBlocked() [3/3]

G4bool G4HadronicInteraction::IsBlocked ( const G4Material aMaterial) const
inherited

Definition at line 193 of file G4HadronicInteraction.cc.

194{
195 for (auto const& mat : theBlockedList) {
196 if (aMaterial == mat) return true;
197 }
198 return false;
199}

References G4HadronicInteraction::theBlockedList.

◆ LowestEnergyLimit()

G4double G4HadronElastic::LowestEnergyLimit ( ) const
inlineinherited

Definition at line 96 of file G4HadronElastic.hh.

97{
98 return lowestEnergyLimit;
99}
G4double lowestEnergyLimit

References G4HadronElastic::lowestEnergyLimit.

Referenced by G4NeutrinoElectronNcModel::ApplyYourself(), and ApplyYourself().

◆ ModelDescription()

void G4NeutronElectronElModel::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronElastic.

Definition at line 101 of file G4NeutronElectronElModel.cc.

102{
103 outFile << "G4NeutronElectronElModel is a neutrino-electron (neutral current) elastic scattering\n"
104 << "model which uses the standard model \n"
105 << "transfer parameterization. The model is fully relativistic\n";
106}

◆ operator!=()

G4bool G4HadronicInteraction::operator!= ( const G4HadronicInteraction right) const
deleteinherited

◆ operator==()

G4bool G4HadronicInteraction::operator== ( const G4HadronicInteraction right) const
deleteinherited

◆ SampleInvariantT()

G4double G4HadronElastic::SampleInvariantT ( const G4ParticleDefinition p,
G4double  plab,
G4int  Z,
G4int  A 
)
overridevirtualinherited

Reimplemented from G4HadronicInteraction.

Reimplemented in G4NuclNuclDiffuseElastic, G4LEHadronProtonElastic, G4LEnp, G4LEpp, G4LowEHadronElastic, and G4hhElastic.

Definition at line 208 of file G4HadronElastic.cc.

210{
211 const G4double plabLowLimit = 400.0*CLHEP::MeV;
212 const G4double GeV2 = GeV*GeV;
213 const G4double z07in13 = std::pow(0.7, 0.3333333333);
214 const G4double numLimit = 18.;
215
216 G4int pdg = std::abs(part->GetPDGEncoding());
217 G4double tmax = pLocalTmax/GeV2;
218
219 G4double aa, bb, cc, dd;
220 G4Pow* g4pow = G4Pow::GetInstance();
221 if (A <= 62) {
222 if (pdg == 211){ //Pions
223 if(mom >= plabLowLimit){ //High energy
224 bb = 14.5*g4pow->Z23(A);/*14.5*/
225 dd = 10.;
226 cc = 0.075*g4pow->Z13(A)/dd;//1.4
227 //aa = g4pow->powZ(A, 1.93)/bb;//1.63
228 aa = (A*A)/bb;//1.63
229 } else { //Low energy
230 bb = 29.*z07in13*z07in13*g4pow->Z23(A);
231 dd = 15.;
232 cc = 0.04*g4pow->Z13(A)/dd;//1.4
233 aa = g4pow->powZ(A, 1.63)/bb;//1.63
234 }
235 } else { //Other particles
236 bb = 14.5*g4pow->Z23(A);
237 dd = 20.;
238 aa = (A*A)/bb;//1.63
239 cc = 1.4*g4pow->Z13(A)/dd;
240 }
241 //===========================
242 } else { //(A>62)
243 if (pdg == 211) {
244 if(mom >= plabLowLimit){ //high
245 bb = 60.*z07in13*g4pow->Z13(A);//60
246 dd = 30.;
247 aa = 0.5*(A*A)/bb;//1.33
248 cc = 4.*g4pow->powZ(A,0.4)/dd;//1:0.4 --- 2: 0.4
249 } else { //low
250 bb = 120.*z07in13*g4pow->Z13(A);//60
251 dd = 30.;
252 aa = 2.*g4pow->powZ(A,1.33)/bb;
253 cc = 4.*g4pow->powZ(A,0.4)/dd;//1:0.4 --- 2: 0.4
254 }
255 } else {
256 bb = 60.*g4pow->Z13(A);
257 dd = 25.;
258 aa = g4pow->powZ(A,1.33)/bb;//1.33
259 cc = 0.2*g4pow->powZ(A,0.4)/dd;//1:0.4 --- 2: 0.4
260 }
261 }
262 G4double q1 = 1.0 - G4Exp(-std::min(bb*tmax, numLimit));
263 G4double q2 = 1.0 - G4Exp(-std::min(dd*tmax, numLimit));
264 G4double s1 = q1*aa;
265 G4double s2 = q2*cc;
266 if((s1 + s2)*G4UniformRand() < s2) {
267 q1 = q2;
268 bb = dd;
269 }
270 return -GeV2*G4Log(1.0 - G4UniformRand()*q1)/bb;
271}
const G4double GeV2
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179
G4double G4Log(G4double x)
Definition: G4Log.hh:226
Definition: G4Pow.hh:49
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double powZ(G4int Z, G4double y) const
Definition: G4Pow.hh:225
G4double Z13(G4int Z) const
Definition: G4Pow.hh:123
G4double Z23(G4int Z) const
Definition: G4Pow.hh:125
static constexpr double MeV
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

References A, G4Exp(), G4Log(), G4UniformRand, G4Pow::GetInstance(), G4ParticleDefinition::GetPDGEncoding(), GeV, GeV2, CLHEP::MeV, G4INCL::Math::min(), G4HadronElastic::pLocalTmax, G4Pow::powZ(), G4Pow::Z13(), and G4Pow::Z23().

Referenced by G4HadronElastic::ApplyYourself(), G4ChipsElasticModel::SampleInvariantT(), G4ElasticHadrNucleusHE::SampleInvariantT(), and G4LowEHadronElastic::SampleInvariantT().

◆ SampleSin2HalfTheta()

G4double G4NeutronElectronElModel::SampleSin2HalfTheta ( G4double  Tkin)

Definition at line 159 of file G4NeutronElectronElModel.cc.

160{
161 G4double result = 0., position;
162 G4int iTkin, iTransfer;
163
164 for( iTkin = 0; iTkin < fEnergyBin; iTkin++)
165 {
166 if( Tkin < fEnergyVector->Energy(iTkin) ) break;
167 }
168 if ( iTkin >= fEnergyBin ) iTkin = fEnergyBin-1; // Tkin is more then theMaxEnergy
169 if ( iTkin < 0 ) iTkin = 0; // against negative index, Tkin < theMinEnergy
170
171 position = (*(*fAngleTable)(iTkin))(fAngleBin-1)*G4UniformRand();
172
173 // G4cout<<"position = "<<position<<G4endl;
174
175 for( iTransfer = 0; iTransfer < fAngleBin; iTransfer++)
176 {
177 if( position <= (*(*fAngleTable)(iTkin))(iTransfer) ) break;
178 }
179 if (iTransfer >= fAngleBin-1) iTransfer = fAngleBin-1;
180
181 // G4cout<<"iTransfer = "<<iTransfer<<G4endl;
182
183 result = GetTransfer(iTkin, iTransfer, position);
184
185 // G4cout<<"t = "<<t<<G4endl;
186
187
188 return result;
189}
G4double GetTransfer(G4int iTkin, G4int iTransfer, G4double position)
#define position
Definition: xmlparse.cc:622

References fAngleBin, fAngleTable, fEnergyBin, G4UniformRand, GetTransfer(), and position.

Referenced by ApplyYourself().

◆ SetCutEnergy()

void G4NeutronElectronElModel::SetCutEnergy ( G4double  ec)
inline

Definition at line 79 of file G4NeutronElectronElModel.hh.

79{fCutEnergy=ec;};

References fCutEnergy.

◆ SetEnergyMomentumCheckLevels()

void G4HadronicInteraction::SetEnergyMomentumCheckLevels ( G4double  relativeLevel,
G4double  absoluteLevel 
)
inlineinherited

Definition at line 149 of file G4HadronicInteraction.hh.

150 { epCheckLevels.first = relativeLevel;
151 epCheckLevels.second = absoluteLevel; }

References G4HadronicInteraction::epCheckLevels.

Referenced by G4BinaryCascade::G4BinaryCascade(), G4CascadeInterface::G4CascadeInterface(), and G4FTFModel::G4FTFModel().

◆ SetLowestEnergyLimit()

void G4HadronElastic::SetLowestEnergyLimit ( G4double  value)
inlineinherited

◆ SetMaxEnergy() [1/3]

void G4HadronicInteraction::SetMaxEnergy ( const G4double  anEnergy)
inlineinherited

Definition at line 102 of file G4HadronicInteraction.hh.

103 { theMaxEnergy = anEnergy; }

References G4HadronicInteraction::theMaxEnergy.

Referenced by G4HadronicInteraction::ActivateFor(), G4IonINCLXXPhysics::AddProcess(), G4BertiniElectroNuclearBuilder::Build(), G4LENDBertiniGammaElectroNuclearBuilder::Build(), G4NeutronLENDBuilder::Build(), G4NeutronPHPBuilder::Build(), G4AlphaPHPBuilder::Build(), G4BertiniKaonBuilder::Build(), G4BertiniNeutronBuilder::Build(), G4BertiniPiKBuilder::Build(), G4BertiniPionBuilder::Build(), G4BertiniProtonBuilder::Build(), G4BinaryAlphaBuilder::Build(), G4BinaryDeuteronBuilder::Build(), G4BinaryHe3Builder::Build(), G4BinaryNeutronBuilder::Build(), G4BinaryPiKBuilder::Build(), G4BinaryPionBuilder::Build(), G4BinaryProtonBuilder::Build(), G4BinaryTritonBuilder::Build(), G4DeuteronPHPBuilder::Build(), G4FTFBinaryKaonBuilder::Build(), G4FTFBinaryNeutronBuilder::Build(), G4FTFBinaryPionBuilder::Build(), G4FTFBinaryProtonBuilder::Build(), G4FTFPAntiBarionBuilder::Build(), G4FTFPKaonBuilder::Build(), G4FTFPNeutronBuilder::Build(), G4FTFPPiKBuilder::Build(), G4FTFPPionBuilder::Build(), G4FTFPProtonBuilder::Build(), G4He3PHPBuilder::Build(), G4HyperonFTFPBuilder::Build(), G4HyperonQGSPBuilder::Build(), G4INCLXXNeutronBuilder::Build(), G4INCLXXPionBuilder::Build(), G4INCLXXProtonBuilder::Build(), G4PrecoNeutronBuilder::Build(), G4PrecoProtonBuilder::Build(), G4ProtonPHPBuilder::Build(), G4QGSBinaryKaonBuilder::Build(), G4QGSBinaryNeutronBuilder::Build(), G4QGSBinaryPiKBuilder::Build(), G4QGSBinaryPionBuilder::Build(), G4QGSBinaryProtonBuilder::Build(), G4QGSPAntiBarionBuilder::Build(), G4QGSPKaonBuilder::Build(), G4QGSPLundStrFragmProtonBuilder::Build(), G4QGSPNeutronBuilder::Build(), G4QGSPPiKBuilder::Build(), G4QGSPPionBuilder::Build(), G4TritonPHPBuilder::Build(), G4QGSPProtonBuilder::Build(), G4HadronicBuilder::BuildFTFP_BERT(), G4HadronicBuilder::BuildFTFQGSP_BERT(), G4QGSBuilder::BuildModel(), G4VHadronPhysics::BuildModel(), G4HadronicBuilder::BuildQGSP_FTFP_BERT(), G4EmExtraPhysics::ConstructGammaElectroNuclear(), LBE::ConstructHad(), G4EmExtraPhysics::ConstructLENDGammaNuclear(), G4HadronDElasticPhysics::ConstructProcess(), G4HadronElasticPhysics::ConstructProcess(), G4HadronHElasticPhysics::ConstructProcess(), G4IonINCLXXPhysics::ConstructProcess(), G4IonPhysics::ConstructProcess(), G4IonPhysicsPHP::ConstructProcess(), G4IonQMDPhysics::ConstructProcess(), G4ANuElNucleusNcModel::G4ANuElNucleusNcModel(), G4ANuMuNucleusNcModel::G4ANuMuNucleusNcModel(), G4BertiniKaonBuilder::G4BertiniKaonBuilder(), G4BertiniPiKBuilder::G4BertiniPiKBuilder(), G4BertiniPionBuilder::G4BertiniPionBuilder(), G4BinaryCascade::G4BinaryCascade(), G4BinaryPiKBuilder::G4BinaryPiKBuilder(), G4BinaryPionBuilder::G4BinaryPionBuilder(), G4ChargeExchange::G4ChargeExchange(), G4DiffuseElastic::G4DiffuseElastic(), G4DiffuseElasticV2::G4DiffuseElasticV2(), G4ElectroVDNuclearModel::G4ElectroVDNuclearModel(), G4EMDissociation::G4EMDissociation(), G4FissLib::G4FissLib(), G4FTFBinaryKaonBuilder::G4FTFBinaryKaonBuilder(), G4FTFBinaryNeutronBuilder::G4FTFBinaryNeutronBuilder(), G4FTFBinaryPiKBuilder::G4FTFBinaryPiKBuilder(), G4FTFBinaryPionBuilder::G4FTFBinaryPionBuilder(), G4FTFBinaryProtonBuilder::G4FTFBinaryProtonBuilder(), G4FTFPAntiBarionBuilder::G4FTFPAntiBarionBuilder(), G4FTFPKaonBuilder::G4FTFPKaonBuilder(), G4FTFPNeutronBuilder::G4FTFPNeutronBuilder(), G4FTFPPiKBuilder::G4FTFPPiKBuilder(), G4FTFPPionBuilder::G4FTFPPionBuilder(), G4FTFPProtonBuilder::G4FTFPProtonBuilder(), G4HadronElastic::G4HadronElastic(), G4HadronicAbsorptionFritiof::G4HadronicAbsorptionFritiof(), G4HadronicAbsorptionFritiofWithBinaryCascade::G4HadronicAbsorptionFritiofWithBinaryCascade(), G4hhElastic::G4hhElastic(), G4HyperonFTFPBuilder::G4HyperonFTFPBuilder(), G4HyperonQGSPBuilder::G4HyperonQGSPBuilder(), G4INCLXXPionBuilder::G4INCLXXPionBuilder(), G4LEHadronProtonElastic::G4LEHadronProtonElastic(), G4LENDModel::G4LENDModel(), G4LEnp::G4LEnp(), G4LEpp::G4LEpp(), G4LFission::G4LFission(), G4LowEGammaNuclearModel::G4LowEGammaNuclearModel(), G4MuonVDNuclearModel::G4MuonVDNuclearModel(), G4NeutrinoElectronCcModel::G4NeutrinoElectronCcModel(), G4NeutrinoElectronNcModel::G4NeutrinoElectronNcModel(), G4NeutrinoNucleusModel::G4NeutrinoNucleusModel(), G4NeutronElectronElModel(), G4NeutronRadCapture::G4NeutronRadCapture(), G4NuclNuclDiffuseElastic::G4NuclNuclDiffuseElastic(), G4NuElNucleusNcModel::G4NuElNucleusNcModel(), G4NuMuNucleusNcModel::G4NuMuNucleusNcModel(), G4ParticleHPCapture::G4ParticleHPCapture(), G4ParticleHPElastic::G4ParticleHPElastic(), G4ParticleHPFission::G4ParticleHPFission(), G4ParticleHPInelastic::G4ParticleHPInelastic(), G4ParticleHPThermalScattering::G4ParticleHPThermalScattering(), G4QGSPAntiBarionBuilder::G4QGSPAntiBarionBuilder(), G4WilsonAbrasionModel::G4WilsonAbrasionModel(), G4HadronPhysicsFTFP_BERT_HP::Neutron(), G4HadronPhysicsINCLXX::Neutron(), G4HadronPhysicsQGSP_BERT_HP::Neutron(), G4HadronPhysicsQGSP_BIC_HP::Neutron(), G4HadronPhysicsShielding::Neutron(), and G4VHadronPhysics::NewModel().

◆ SetMaxEnergy() [2/3]

void G4HadronicInteraction::SetMaxEnergy ( G4double  anEnergy,
const G4Element anElement 
)
inherited

Definition at line 151 of file G4HadronicInteraction.cc.

153{
154 Block();
155 if(!theMaxEnergyListElements.empty()) {
156 for(auto & elmlist : theMaxEnergyListElements) {
157 if( anElement == elmlist.second ) {
158 elmlist.first = anEnergy;
159 return;
160 }
161 }
162 }
163 theMaxEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
164}

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theMaxEnergyListElements.

◆ SetMaxEnergy() [3/3]

void G4HadronicInteraction::SetMaxEnergy ( G4double  anEnergy,
const G4Material aMaterial 
)
inherited

Definition at line 166 of file G4HadronicInteraction.cc.

167{
168 Block();
169 if(!theMaxEnergyList.empty()) {
170 for(auto & matlist: theMaxEnergyList) {
171 if( aMaterial == matlist.second ) {
172 matlist.first = anEnergy;
173 return;
174 }
175 }
176 }
177 theMaxEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
178}

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theMaxEnergyList.

◆ SetMinEnergy() [1/3]

void G4HadronicInteraction::SetMinEnergy ( G4double  anEnergy)
inlineinherited

Definition at line 89 of file G4HadronicInteraction.hh.

90 { theMinEnergy = anEnergy; }

References G4HadronicInteraction::theMinEnergy.

Referenced by G4HadronicInteraction::ActivateFor(), G4BertiniElectroNuclearBuilder::Build(), G4LENDBertiniGammaElectroNuclearBuilder::Build(), G4NeutronLENDBuilder::Build(), G4NeutronPHPBuilder::Build(), G4AlphaPHPBuilder::Build(), G4BertiniKaonBuilder::Build(), G4BertiniNeutronBuilder::Build(), G4BertiniPiKBuilder::Build(), G4BertiniPionBuilder::Build(), G4BertiniProtonBuilder::Build(), G4BinaryAlphaBuilder::Build(), G4BinaryDeuteronBuilder::Build(), G4BinaryHe3Builder::Build(), G4BinaryNeutronBuilder::Build(), G4BinaryPiKBuilder::Build(), G4BinaryPionBuilder::Build(), G4BinaryProtonBuilder::Build(), G4BinaryTritonBuilder::Build(), G4DeuteronPHPBuilder::Build(), G4FTFBinaryKaonBuilder::Build(), G4FTFBinaryNeutronBuilder::Build(), G4FTFBinaryPiKBuilder::Build(), G4FTFBinaryPionBuilder::Build(), G4FTFBinaryProtonBuilder::Build(), G4FTFPAntiBarionBuilder::Build(), G4FTFPKaonBuilder::Build(), G4FTFPNeutronBuilder::Build(), G4FTFPPiKBuilder::Build(), G4FTFPPionBuilder::Build(), G4FTFPProtonBuilder::Build(), G4He3PHPBuilder::Build(), G4HyperonFTFPBuilder::Build(), G4HyperonQGSPBuilder::Build(), G4INCLXXNeutronBuilder::Build(), G4INCLXXPionBuilder::Build(), G4INCLXXProtonBuilder::Build(), G4PrecoNeutronBuilder::Build(), G4PrecoProtonBuilder::Build(), G4ProtonPHPBuilder::Build(), G4QGSBinaryKaonBuilder::Build(), G4QGSBinaryNeutronBuilder::Build(), G4QGSBinaryPiKBuilder::Build(), G4QGSBinaryPionBuilder::Build(), G4QGSBinaryProtonBuilder::Build(), G4QGSPAntiBarionBuilder::Build(), G4QGSPKaonBuilder::Build(), G4QGSPLundStrFragmProtonBuilder::Build(), G4QGSPNeutronBuilder::Build(), G4QGSPPiKBuilder::Build(), G4QGSPPionBuilder::Build(), G4TritonPHPBuilder::Build(), G4QGSPProtonBuilder::Build(), G4QGSBuilder::BuildModel(), G4VHadronPhysics::BuildModel(), G4EmExtraPhysics::ConstructGammaElectroNuclear(), LBE::ConstructHad(), G4EmExtraPhysics::ConstructLENDGammaNuclear(), G4HadronElasticPhysicsHP::ConstructProcess(), G4HadronElasticPhysicsLEND::ConstructProcess(), G4HadronElasticPhysicsPHP::ConstructProcess(), G4HadronDElasticPhysics::ConstructProcess(), G4HadronElasticPhysics::ConstructProcess(), G4HadronHElasticPhysics::ConstructProcess(), G4IonElasticPhysics::ConstructProcess(), G4IonINCLXXPhysics::ConstructProcess(), G4IonPhysics::ConstructProcess(), G4IonPhysicsPHP::ConstructProcess(), G4IonQMDPhysics::ConstructProcess(), G4ANuElNucleusNcModel::G4ANuElNucleusNcModel(), G4ANuMuNucleusNcModel::G4ANuMuNucleusNcModel(), G4BertiniKaonBuilder::G4BertiniKaonBuilder(), G4BertiniPiKBuilder::G4BertiniPiKBuilder(), G4BertiniPionBuilder::G4BertiniPionBuilder(), G4BinaryCascade::G4BinaryCascade(), G4BinaryPiKBuilder::G4BinaryPiKBuilder(), G4BinaryPionBuilder::G4BinaryPionBuilder(), G4ChargeExchange::G4ChargeExchange(), G4DiffuseElastic::G4DiffuseElastic(), G4DiffuseElasticV2::G4DiffuseElasticV2(), G4ElectroVDNuclearModel::G4ElectroVDNuclearModel(), G4EMDissociation::G4EMDissociation(), G4FissLib::G4FissLib(), G4FTFBinaryKaonBuilder::G4FTFBinaryKaonBuilder(), G4FTFBinaryNeutronBuilder::G4FTFBinaryNeutronBuilder(), G4FTFBinaryPiKBuilder::G4FTFBinaryPiKBuilder(), G4FTFBinaryPionBuilder::G4FTFBinaryPionBuilder(), G4FTFBinaryProtonBuilder::G4FTFBinaryProtonBuilder(), G4FTFPAntiBarionBuilder::G4FTFPAntiBarionBuilder(), G4FTFPKaonBuilder::G4FTFPKaonBuilder(), G4FTFPNeutronBuilder::G4FTFPNeutronBuilder(), G4FTFPPiKBuilder::G4FTFPPiKBuilder(), G4FTFPPionBuilder::G4FTFPPionBuilder(), G4FTFPProtonBuilder::G4FTFPProtonBuilder(), G4HadronElastic::G4HadronElastic(), G4HadronicAbsorptionBertini::G4HadronicAbsorptionBertini(), G4HadronicAbsorptionFritiof::G4HadronicAbsorptionFritiof(), G4HadronicAbsorptionFritiofWithBinaryCascade::G4HadronicAbsorptionFritiofWithBinaryCascade(), G4hhElastic::G4hhElastic(), G4HyperonFTFPBuilder::G4HyperonFTFPBuilder(), G4HyperonQGSPBuilder::G4HyperonQGSPBuilder(), G4INCLXXPionBuilder::G4INCLXXPionBuilder(), G4LEHadronProtonElastic::G4LEHadronProtonElastic(), G4LENDModel::G4LENDModel(), G4LEnp::G4LEnp(), G4LEpp::G4LEpp(), G4LFission::G4LFission(), G4LowEGammaNuclearModel::G4LowEGammaNuclearModel(), G4MuonVDNuclearModel::G4MuonVDNuclearModel(), G4NeutrinoElectronCcModel::G4NeutrinoElectronCcModel(), G4NeutrinoElectronNcModel::G4NeutrinoElectronNcModel(), G4NeutrinoNucleusModel::G4NeutrinoNucleusModel(), G4NeutronElectronElModel(), G4NeutronRadCapture::G4NeutronRadCapture(), G4NuclNuclDiffuseElastic::G4NuclNuclDiffuseElastic(), G4NuElNucleusNcModel::G4NuElNucleusNcModel(), G4NuMuNucleusNcModel::G4NuMuNucleusNcModel(), G4ParticleHPCapture::G4ParticleHPCapture(), G4ParticleHPElastic::G4ParticleHPElastic(), G4ParticleHPFission::G4ParticleHPFission(), G4ParticleHPInelastic::G4ParticleHPInelastic(), G4ParticleHPThermalScattering::G4ParticleHPThermalScattering(), G4QGSPAntiBarionBuilder::G4QGSPAntiBarionBuilder(), G4WilsonAbrasionModel::G4WilsonAbrasionModel(), G4NeutrinoElectronCcModel::IsApplicable(), G4HadronPhysicsFTFP_BERT_HP::Neutron(), G4HadronPhysicsINCLXX::Neutron(), G4HadronPhysicsQGSP_BERT_HP::Neutron(), G4HadronPhysicsQGSP_BIC_HP::Neutron(), G4HadronPhysicsShielding::Neutron(), and G4VHadronPhysics::NewModel().

◆ SetMinEnergy() [2/3]

void G4HadronicInteraction::SetMinEnergy ( G4double  anEnergy,
const G4Element anElement 
)
inherited

Definition at line 101 of file G4HadronicInteraction.cc.

103{
104 Block();
105 if(!theMinEnergyListElements.empty()) {
106 for(auto & elmlist : theMinEnergyListElements) {
107 if( anElement == elmlist.second ) {
108 elmlist.first = anEnergy;
109 return;
110 }
111 }
112 }
113 theMinEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
114}

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theMinEnergyListElements.

◆ SetMinEnergy() [3/3]

void G4HadronicInteraction::SetMinEnergy ( G4double  anEnergy,
const G4Material aMaterial 
)
inherited

Definition at line 116 of file G4HadronicInteraction.cc.

118{
119 Block();
120 if(!theMinEnergyList.empty()) {
121 for(auto & matlist : theMinEnergyList) {
122 if( aMaterial == matlist.second ) {
123 matlist.first = anEnergy;
124 return;
125 }
126 }
127 }
128 theMinEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
129}

References G4HadronicInteraction::Block(), and G4HadronicInteraction::theMinEnergyList.

◆ SetModelName()

void G4HadronicInteraction::SetModelName ( const G4String nam)
inlineprotectedinherited

◆ SetRecoilEnergyThreshold()

void G4HadronicInteraction::SetRecoilEnergyThreshold ( G4double  val)
inlineinherited

◆ SetVerboseLevel()

void G4HadronicInteraction::SetVerboseLevel ( G4int  value)
inlineinherited

◆ XscIntegrand()

G4double G4NeutronElectronElModel::XscIntegrand ( G4double  x)

Definition at line 239 of file G4NeutronElectronElModel.cc.

240{
241 G4double result = 1., q2, znq2, znf, znf2, znf4;
242
243 znq2 = 1. + 2.*fee*x/fM;
244
245 q2 = 4.*fee2*x/znq2;
246
247 znf = 1 + q2/fMv2;
248 znf2 = znf*znf;
249 znf4 = znf2*znf2;
250
251 result /= ( x + fAm )*znq2*znq2*znf4;
252
253 result *= ( 1 - x )/( 1 + q2/4./fM2 ) + 2.*x;
254
255 return result;
256}

References fAm, fee, fee2, fM, fM2, and fMv2.

Referenced by Initialise().

Field Documentation

◆ epCheckLevels

std::pair<G4double, G4double> G4HadronicInteraction::epCheckLevels
privateinherited

◆ fAm

G4double G4NeutronElectronElModel::fAm
private

◆ fAngleBin

G4int G4NeutronElectronElModel::fAngleBin
private

◆ fAngleTable

G4PhysicsTable* G4NeutronElectronElModel::fAngleTable
private

◆ fCutEnergy

G4double G4NeutronElectronElModel::fCutEnergy
private

◆ fee

G4double G4NeutronElectronElModel::fee
private

Definition at line 88 of file G4NeutronElectronElModel.hh.

Referenced by ApplyYourself(), CalculateAm(), and XscIntegrand().

◆ fee2

G4double G4NeutronElectronElModel::fee2
private

Definition at line 88 of file G4NeutronElectronElModel.hh.

Referenced by CalculateAm(), and XscIntegrand().

◆ fEnergyBin

G4int G4NeutronElectronElModel::fEnergyBin
private

◆ fEnergyVector

G4PhysicsLogVector* G4NeutronElectronElModel::fEnergyVector
private

◆ fM

G4double G4NeutronElectronElModel::fM
private

◆ fM2

G4double G4NeutronElectronElModel::fM2
private

Definition at line 88 of file G4NeutronElectronElModel.hh.

Referenced by G4NeutronElectronElModel(), and XscIntegrand().

◆ fMaxEnergy

G4double G4NeutronElectronElModel::fMaxEnergy
private

Definition at line 89 of file G4NeutronElectronElModel.hh.

Referenced by G4NeutronElectronElModel(), and IsApplicable().

◆ fme

G4double G4NeutronElectronElModel::fme
private

◆ fme2

G4double G4NeutronElectronElModel::fme2
private

Definition at line 88 of file G4NeutronElectronElModel.hh.

Referenced by CalculateAm(), and G4NeutronElectronElModel().

◆ fMinEnergy

G4double G4NeutronElectronElModel::fMinEnergy
private

Definition at line 89 of file G4NeutronElectronElModel.hh.

Referenced by G4NeutronElectronElModel(), and IsApplicable().

◆ fMv2

G4double G4NeutronElectronElModel::fMv2
private

Definition at line 88 of file G4NeutronElectronElModel.hh.

Referenced by G4NeutronElectronElModel(), and XscIntegrand().

◆ isBlocked

G4bool G4HadronicInteraction::isBlocked
protectedinherited

◆ lowestEnergyLimit

G4double G4HadronElastic::lowestEnergyLimit
privateinherited

◆ nwarn

G4int G4HadronElastic::nwarn
privateinherited

◆ pLocalTmax

G4double G4HadronElastic::pLocalTmax
protectedinherited

◆ recoilEnergyThreshold

G4double G4HadronicInteraction::recoilEnergyThreshold
privateinherited

◆ registry

G4HadronicInteractionRegistry* G4HadronicInteraction::registry
privateinherited

◆ secID

G4int G4HadronElastic::secID
protectedinherited

◆ theAlpha

G4ParticleDefinition* G4HadronElastic::theAlpha
privateinherited

◆ theBlockedList

std::vector<const G4Material *> G4HadronicInteraction::theBlockedList
privateinherited

◆ theBlockedListElements

std::vector<const G4Element *> G4HadronicInteraction::theBlockedListElements
privateinherited

◆ theDeuteron

G4ParticleDefinition* G4HadronElastic::theDeuteron
privateinherited

◆ theElectron

G4ParticleDefinition* G4NeutronElectronElModel::theElectron
private

Definition at line 91 of file G4NeutronElectronElModel.hh.

Referenced by ApplyYourself(), and G4NeutronElectronElModel().

◆ theMaxEnergy

G4double G4HadronicInteraction::theMaxEnergy
protectedinherited

◆ theMaxEnergyList

std::vector<std::pair<G4double, const G4Material *> > G4HadronicInteraction::theMaxEnergyList
privateinherited

◆ theMaxEnergyListElements

std::vector<std::pair<G4double, const G4Element *> > G4HadronicInteraction::theMaxEnergyListElements
privateinherited

◆ theMinEnergy

G4double G4HadronicInteraction::theMinEnergy
protectedinherited

◆ theMinEnergyList

std::vector<std::pair<G4double, const G4Material *> > G4HadronicInteraction::theMinEnergyList
privateinherited

◆ theMinEnergyListElements

std::vector<std::pair<G4double, const G4Element *> > G4HadronicInteraction::theMinEnergyListElements
privateinherited

◆ theModelName

G4String G4HadronicInteraction::theModelName
privateinherited

◆ theNeutron

G4ParticleDefinition* G4HadronElastic::theNeutron
privateinherited

Definition at line 83 of file G4HadronElastic.hh.

Referenced by G4HadronElastic::G4HadronElastic().

◆ theParticleChange

G4HadFinalState G4HadronicInteraction::theParticleChange
protectedinherited

Definition at line 172 of file G4HadronicInteraction.hh.

Referenced by G4WilsonAbrasionModel::ApplyYourself(), G4EMDissociation::ApplyYourself(), G4LENDCapture::ApplyYourself(), G4LENDElastic::ApplyYourself(), G4LENDFission::ApplyYourself(), G4LENDInelastic::ApplyYourself(), G4ElectroVDNuclearModel::ApplyYourself(), G4ParticleHPThermalScattering::ApplyYourself(), G4NeutrinoElectronNcModel::ApplyYourself(), ApplyYourself(), G4LFission::ApplyYourself(), G4ANuElNucleusCcModel::ApplyYourself(), G4ANuElNucleusNcModel::ApplyYourself(), G4ANuMuNucleusCcModel::ApplyYourself(), G4ANuMuNucleusNcModel::ApplyYourself(), G4MuonVDNuclearModel::ApplyYourself(), G4NeutrinoElectronCcModel::ApplyYourself(), G4NuElNucleusCcModel::ApplyYourself(), G4NuElNucleusNcModel::ApplyYourself(), G4NuMuNucleusCcModel::ApplyYourself(), G4NuMuNucleusNcModel::ApplyYourself(), G4QMDReaction::ApplyYourself(), G4NeutronRadCapture::ApplyYourself(), G4LowEGammaNuclearModel::ApplyYourself(), G4ChargeExchange::ApplyYourself(), G4HadronElastic::ApplyYourself(), G4LEHadronProtonElastic::ApplyYourself(), G4LEnp::ApplyYourself(), G4LEpp::ApplyYourself(), G4BinaryCascade::ApplyYourself(), G4CascadeInterface::ApplyYourself(), G4LMsdGenerator::ApplyYourself(), G4ElectroVDNuclearModel::CalculateEMVertex(), G4MuonVDNuclearModel::CalculateEMVertex(), G4ElectroVDNuclearModel::CalculateHadronicVertex(), G4MuonVDNuclearModel::CalculateHadronicVertex(), G4NeutrinoNucleusModel::CoherentPion(), G4CascadeInterface::copyOutputToHadronicResult(), G4BinaryCascade::DebugEpConservation(), G4BinaryCascade::DebugFinalEpConservation(), G4NeutrinoNucleusModel::FinalBarion(), G4NeutrinoNucleusModel::FinalMeson(), G4WilsonAbrasionModel::GetAbradedNucleons(), G4CascadeInterface::NoInteraction(), G4CascadeInterface::Propagate(), G4NeutrinoNucleusModel::RecoilDeexcitation(), G4LEHadronProtonElastic::~G4LEHadronProtonElastic(), G4LEnp::~G4LEnp(), and G4LFission::~G4LFission().

◆ theProton

G4ParticleDefinition* G4HadronElastic::theProton
privateinherited

◆ verboseLevel

G4int G4HadronicInteraction::verboseLevel
protectedinherited

Definition at line 177 of file G4HadronicInteraction.hh.

Referenced by G4WilsonAbrasionModel::ApplyYourself(), G4EMDissociation::ApplyYourself(), G4LFission::ApplyYourself(), G4MuMinusCapturePrecompound::ApplyYourself(), G4NeutronRadCapture::ApplyYourself(), G4LowEGammaNuclearModel::ApplyYourself(), G4ChargeExchange::ApplyYourself(), G4HadronElastic::ApplyYourself(), G4LEHadronProtonElastic::ApplyYourself(), G4LEnp::ApplyYourself(), G4LEpp::ApplyYourself(), G4CascadeInterface::ApplyYourself(), G4CascadeInterface::checkFinalResult(), G4CascadeInterface::copyOutputToHadronicResult(), G4CascadeInterface::copyOutputToReactionProducts(), G4LENDModel::create_used_target_map(), G4CascadeInterface::createBullet(), G4CascadeInterface::createTarget(), G4ElasticHadrNucleusHE::DefineHadronValues(), G4ElasticHadrNucleusHE::FillData(), G4ElasticHadrNucleusHE::FillFq2(), G4DiffuseElastic::G4DiffuseElastic(), G4DiffuseElasticV2::G4DiffuseElasticV2(), G4ElasticHadrNucleusHE::G4ElasticHadrNucleusHE(), G4EMDissociation::G4EMDissociation(), G4hhElastic::G4hhElastic(), G4NuclNuclDiffuseElastic::G4NuclNuclDiffuseElastic(), G4WilsonAbrasionModel::G4WilsonAbrasionModel(), G4ElasticHadrNucleusHE::GetFt(), G4ElasticHadrNucleusHE::GetLightFq2(), G4ElasticHadrNucleusHE::GetQ2_2(), G4HadronicInteraction::GetVerboseLevel(), G4ElasticHadrNucleusHE::HadronNucleusQ2_2(), G4ElasticHadrNucleusHE::HadronProtonQ2(), G4LFission::init(), G4DiffuseElastic::Initialise(), G4DiffuseElasticV2::Initialise(), G4NuclNuclDiffuseElastic::Initialise(), G4DiffuseElastic::InitialiseOnFly(), G4DiffuseElasticV2::InitialiseOnFly(), G4NuclNuclDiffuseElastic::InitialiseOnFly(), G4CascadeInterface::makeDynamicParticle(), G4CascadeInterface::NoInteraction(), G4CascadeInterface::Propagate(), G4ElasticHadrNucleusHE::SampleInvariantT(), G4AntiNuclElastic::SampleThetaCMS(), G4DiffuseElastic::SampleThetaLab(), G4NuclNuclDiffuseElastic::SampleThetaLab(), G4AntiNuclElastic::SampleThetaLab(), G4WilsonAbrasionModel::SetUseAblation(), G4HadronicInteraction::SetVerboseLevel(), G4WilsonAbrasionModel::SetVerboseLevel(), G4DiffuseElastic::ThetaCMStoThetaLab(), G4DiffuseElasticV2::ThetaCMStoThetaLab(), G4NuclNuclDiffuseElastic::ThetaCMStoThetaLab(), G4DiffuseElastic::ThetaLabToThetaCMS(), G4DiffuseElasticV2::ThetaLabToThetaCMS(), and G4NuclNuclDiffuseElastic::ThetaLabToThetaCMS().


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