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

#include <G4MuonVDNuclearModel.hh>

Inheritance diagram for G4MuonVDNuclearModel:
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 &)
 
void DeActivateFor (const G4Element *anElement)
 
void DeActivateFor (const G4Material *aMaterial)
 
 G4MuonVDNuclearModel ()
 
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
 
G4int GetVerboseLevel () const
 
virtual void InitialiseModel ()
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4bool IsBlocked (const G4Element *anElement) const
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
virtual void ModelDescription (std::ostream &outFile) const
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
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)
 
virtual ~G4MuonVDNuclearModel ()
 

Protected Member Functions

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

Protected Attributes

G4bool isBlocked
 
G4double theMaxEnergy
 
G4double theMinEnergy
 
G4HadFinalState theParticleChange
 
G4int verboseLevel
 

Private Member Functions

G4DynamicParticleCalculateEMVertex (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
void CalculateHadronicVertex (G4DynamicParticle *incident, G4Nucleus &target)
 
 G4MuonVDNuclearModel (const G4MuonVDNuclearModel &)=delete
 
void MakeSamplingTable ()
 
G4MuonVDNuclearModeloperator= (const G4MuonVDNuclearModel &right)=delete
 

Private Attributes

G4CascadeInterfacebert
 
G4double CutFixed
 
std::pair< G4double, G4doubleepCheckLevels
 
G4TheoFSGeneratorftfp
 
G4bool isMaster
 
G4KokoulinMuonNuclearXSmuNucXS
 
G4double recoilEnergyThreshold
 
G4HadronicInteractionRegistryregistry
 
G4int secID
 
std::vector< const G4Material * > theBlockedList
 
std::vector< const G4Element * > theBlockedListElements
 
G4LundStringFragmentationtheFragmentation
 
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
 
G4ExcitedStringDecaytheStringDecay
 

Static Private Attributes

static const G4double adat [nzdat] = {1.01,9.01,26.98,63.55,238.03}
 
static G4ElementDatafElementData = nullptr
 
static const G4int NBIN = 800
 
static const G4int ntdat = 73
 
static const G4int nzdat = 5
 
static const G4double tdat [ntdat]
 
static const G4int zdat [nzdat] = {1, 4, 13, 29, 92}
 

Detailed Description

Definition at line 49 of file G4MuonVDNuclearModel.hh.

Constructor & Destructor Documentation

◆ G4MuonVDNuclearModel() [1/2]

G4MuonVDNuclearModel::G4MuonVDNuclearModel ( )
explicit

Definition at line 72 of file G4MuonVDNuclearModel.cc.

73 : G4HadronicInteraction("G4MuonVDNuclearModel"),isMaster(false)
74{
76 GetCrossSectionDataSet(G4KokoulinMuonNuclearXS::Default_Name());
77
78 SetMinEnergy(0.0);
80 CutFixed = 0.2*CLHEP::GeV;
81
85 isMaster = true;
86 }
87
88 // reuse existing pre-compound model
89 G4GeneratorPrecompoundInterface* precoInterface
93 G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(p);
94 if(!pre) { pre = new G4PreCompoundModel(); }
95 precoInterface->SetDeExcitation(pre);
96
97 // Build FTFP model
98 ftfp = new G4TheoFSGenerator();
99 ftfp->SetTransport(precoInterface);
102 G4FTFModel* theStringModel = new G4FTFModel;
103 theStringModel->SetFragmentationModel(theStringDecay);
104 ftfp->SetHighEnergyGenerator(theStringModel);
105
106 // Build Bertini cascade
107 bert = new G4CascadeInterface();
108
109 // Creator model ID
111}
static G4CrossSectionDataSetRegistry * Instance()
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
void SetMinEnergy(G4double anEnergy)
G4HadronicInteraction(const G4String &modelName="HadronicModel")
const G4String & GetModelName() const
void SetMaxEnergy(const G4double anEnergy)
static const char * Default_Name()
G4TheoFSGenerator * ftfp
G4ExcitedStringDecay * theStringDecay
G4LundStringFragmentation * theFragmentation
static G4ElementData * fElementData
G4KokoulinMuonNuclearXS * muNucXS
G4CascadeInterface * bert
static G4int GetModelID(const G4int modelIndex)
void SetTransport(G4VIntraNuclearTransportModel *const value)
void SetHighEnergyGenerator(G4VHighEnergyGenerator *const value)
void SetDeExcitation(G4VPreCompoundModel *ptr)
void SetFragmentationModel(G4VStringFragmentation *aModel)
static constexpr double GeV
static constexpr double PeV
G4bool IsMasterThread()
Definition: G4Threading.cc:124

References bert, CutFixed, G4KokoulinMuonNuclearXS::Default_Name(), fElementData, G4HadronicInteractionRegistry::FindModel(), ftfp, G4PhysicsModelCatalog::GetModelID(), G4HadronicInteraction::GetModelName(), CLHEP::GeV, G4CrossSectionDataSetRegistry::Instance(), G4HadronicInteractionRegistry::Instance(), isMaster, G4Threading::IsMasterThread(), MakeSamplingTable(), muNucXS, CLHEP::PeV, secID, G4VIntraNuclearTransportModel::SetDeExcitation(), G4VPartonStringModel::SetFragmentationModel(), G4TheoFSGenerator::SetHighEnergyGenerator(), G4HadronicInteraction::SetMaxEnergy(), G4HadronicInteraction::SetMinEnergy(), G4TheoFSGenerator::SetTransport(), theFragmentation, and theStringDecay.

◆ ~G4MuonVDNuclearModel()

G4MuonVDNuclearModel::~G4MuonVDNuclearModel ( )
virtual

Definition at line 113 of file G4MuonVDNuclearModel.cc.

114{
115 delete theFragmentation;
116 delete theStringDecay;
117
118 if(isMaster) {
119 delete fElementData;
120 fElementData = nullptr;
121 }
122}

References fElementData, isMaster, theFragmentation, and theStringDecay.

◆ G4MuonVDNuclearModel() [2/2]

G4MuonVDNuclearModel::G4MuonVDNuclearModel ( const G4MuonVDNuclearModel )
privatedelete

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 * G4MuonVDNuclearModel::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 125 of file G4MuonVDNuclearModel.cc.

127{
129
130 // For very low energy, return initial track
131 G4double epmax = aTrack.GetTotalEnergy() - 0.5*proton_mass_c2;
132 if (epmax <= CutFixed) {
136 return &theParticleChange;
137 }
138
139 // Produce recoil muon and transferred photon
140 G4DynamicParticle* transferredPhoton = CalculateEMVertex(aTrack, targetNucleus);
141
142 // Interact the gamma with the nucleus
143 CalculateHadronicVertex(transferredPhoton, targetNucleus);
144 return &theParticleChange;
145}
@ isAlive
double G4double
Definition: G4Types.hh:83
Hep3Vector unit() const
Hep3Vector vect() const
void SetStatusChange(G4HadFinalStateStatus aS)
void SetEnergyChange(G4double anEnergy)
void SetMomentumChange(const G4ThreeVector &aV)
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
G4double GetTotalEnergy() const
void CalculateHadronicVertex(G4DynamicParticle *incident, G4Nucleus &target)
G4DynamicParticle * CalculateEMVertex(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
float proton_mass_c2
Definition: hepunit.py:274

References CalculateEMVertex(), CalculateHadronicVertex(), G4HadFinalState::Clear(), CutFixed, G4HadProjectile::Get4Momentum(), G4HadProjectile::GetKineticEnergy(), G4HadProjectile::GetTotalEnergy(), isAlive, source.hepunit::proton_mass_c2, G4HadFinalState::SetEnergyChange(), G4HadFinalState::SetMomentumChange(), G4HadFinalState::SetStatusChange(), G4HadronicInteraction::theParticleChange, CLHEP::Hep3Vector::unit(), and CLHEP::HepLorentzVector::vect().

◆ Block()

void G4HadronicInteraction::Block ( )
inlineprotectedinherited

◆ BuildPhysicsTable()

void G4HadronicInteraction::BuildPhysicsTable ( const G4ParticleDefinition )
virtualinherited

◆ CalculateEMVertex()

G4DynamicParticle * G4MuonVDNuclearModel::CalculateEMVertex ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
private

Definition at line 148 of file G4MuonVDNuclearModel.cc.

150{
151 // Select sampling table
153 G4double TotalEnergy = aTrack.GetTotalEnergy();
155 G4Pow* g4calc = G4Pow::GetInstance();
156 G4double lnZ = g4calc->logZ(targetNucleus.GetZ_asInt());
157
158 G4double epmin = CutFixed;
159 G4double epmax = TotalEnergy - 0.5*proton_mass_c2;
160 G4double m0 = CutFixed;
161
162 G4double delmin = 1.e10;
163 G4double del;
164 G4int izz = 0;
165 G4int itt = 0;
166
167 for (G4int iz = 0; iz < nzdat; ++iz) {
168 del = std::abs(lnZ - g4calc->logZ(zdat[iz]));
169 if (del < delmin) {
170 delmin = del;
171 izz = iz;
172 }
173 }
174
175 delmin = 1.e10;
176 for (G4int it = 0; it < ntdat; ++it) {
177 del = std::abs(G4Log(KineticEnergy)-G4Log(tdat[it]) );
178 if (del < delmin) {
179 delmin = del;
180 itt = it;
181 }
182 }
183
184 // Sample the energy transfer according to the probability table
186
187 G4int iy;
188
189 G4int Z = zdat[izz];
190
191 for(iy = 0; iy<NBIN; ++iy) {
192
194 if(pvv >= r) { break; }
195 }
196
197 // Sampling is done uniformly in y in the bin
200 G4double y = pvx + G4UniformRand() * (pvx1 - pvx);
201
202 G4double x = G4Exp(y);
203 G4double ep = epmin*G4Exp(x*G4Log(epmax/epmin) );
204
205 // Sample scattering angle of mu, but first t should be sampled.
206 G4double yy = ep/TotalEnergy;
207 G4double tmin = Mass*Mass*yy*yy/(1.-yy);
208 G4double tmax = 2.*proton_mass_c2*ep;
209 G4double t1;
210 G4double t2;
211 if (m0 < ep) {
212 t1 = m0*m0;
213 t2 = ep*ep;
214 } else {
215 t1 = ep*ep;
216 t2 = m0*m0;
217 }
218
219 G4double w1 = tmax*t1;
220 G4double w2 = tmax+t1;
221 G4double w3 = tmax*(tmin+t1)/(tmin*w2);
222 G4double y1 = 1.-yy;
223 G4double y2 = 0.5*yy*yy;
224 G4double y3 = y1+y2;
225
226 G4double t = 0.0;
227 G4double rej = 0.0;
228
229 // Now sample t
230 G4int ntry = 0;
231 do
232 {
233 ntry += 1;
234 if (ntry > 10000) {
236 eda << " While count exceeded " << G4endl;
237 G4Exception("G4MuonVDNuclearModel::CalculateEMVertex()", "HAD_RPG_100", JustWarning, eda);
238 break;
239 }
240
241 t = w1/(w2*G4Exp(G4UniformRand()*G4Log(w3))-tmax);
242 rej = (1.-t/tmax)*(y1*(1.-tmin/t)+y2)/(y3*(1.-t/t2));
243 } while (G4UniformRand() > rej) ; /* Loop checking, 01.09.2015, D.Wright */
244
245 // compute angle from t
246 G4double sinth2 =
247 0.5*(t-tmin)/(2.*(TotalEnergy*(TotalEnergy-ep)-Mass*Mass)-tmin);
248 G4double theta = std::acos(1. - 2.*sinth2);
249
251 G4double sinth = std::sin(theta);
252 G4double dirx = sinth*std::cos(phi);
253 G4double diry = sinth*std::sin(phi);
254 G4double dirz = std::cos(theta);
255 G4ThreeVector finalDirection(dirx,diry,dirz);
256 G4ThreeVector ParticleDirection(aTrack.Get4Momentum().vect().unit() );
257 finalDirection.rotateUz(ParticleDirection);
258
259 G4double NewKinEnergy = KineticEnergy - ep;
260 G4double finalMomentum = std::sqrt(NewKinEnergy*(NewKinEnergy+2.*Mass) );
261 G4double Ef = NewKinEnergy + Mass;
262 G4double initMomentum = std::sqrt(KineticEnergy*(TotalEnergy+Mass) );
263
264 // Set energy and direction of scattered primary in theParticleChange
267 theParticleChange.SetMomentumChange(finalDirection);
268
269 // Now create the emitted gamma
270 G4LorentzVector primaryMomentum(initMomentum*ParticleDirection, TotalEnergy);
271 G4LorentzVector fsMomentum(finalMomentum*finalDirection, Ef);
272 G4LorentzVector momentumTransfer = primaryMomentum - fsMomentum;
273
274 G4DynamicParticle* gamma =
275 new G4DynamicParticle(G4Gamma::Gamma(), momentumTransfer);
276
277 return gamma;
278}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:179
G4double G4Log(G4double x)
Definition: G4Log.hh:226
static constexpr double twopi
Definition: G4SIunits.hh:56
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
#define G4endl
Definition: G4ios.hh:57
#define G4UniformRand()
Definition: Randomize.hh:52
G4Physics2DVector * GetElement2DData(G4int Z)
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
static G4MuonMinus * MuonMinus()
Definition: G4MuonMinus.cc:99
static const G4int NBIN
static const G4int nzdat
static const G4int ntdat
static const G4int zdat[nzdat]
static const G4double tdat[ntdat]
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:105
G4double GetValue(std::size_t idx, std::size_t idy) const
G4double GetX(std::size_t index) const
Definition: G4Pow.hh:49
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double logZ(G4int Z) const
Definition: G4Pow.hh:137

References CutFixed, fElementData, G4endl, G4Exception(), G4Exp(), G4Log(), G4UniformRand, G4Gamma::Gamma(), G4HadProjectile::Get4Momentum(), G4ElementData::GetElement2DData(), G4Pow::GetInstance(), G4HadProjectile::GetKineticEnergy(), G4ParticleDefinition::GetPDGMass(), G4HadProjectile::GetTotalEnergy(), G4Physics2DVector::GetValue(), G4Physics2DVector::GetX(), G4Nucleus::GetZ_asInt(), isAlive, JustWarning, G4Pow::logZ(), G4MuonMinus::MuonMinus(), NBIN, ntdat, nzdat, source.hepunit::proton_mass_c2, CLHEP::Hep3Vector::rotateUz(), G4HadFinalState::SetEnergyChange(), G4HadFinalState::SetMomentumChange(), G4HadFinalState::SetStatusChange(), tdat, G4HadronicInteraction::theParticleChange, twopi, CLHEP::Hep3Vector::unit(), CLHEP::HepLorentzVector::vect(), Z, and zdat.

Referenced by ApplyYourself().

◆ CalculateHadronicVertex()

void G4MuonVDNuclearModel::CalculateHadronicVertex ( G4DynamicParticle incident,
G4Nucleus target 
)
private

Definition at line 281 of file G4MuonVDNuclearModel.cc.

283{
284 G4HadFinalState* hfs = 0;
285 G4double gammaE = incident->GetTotalEnergy();
286
287 if (gammaE < 10*GeV) {
288 G4HadProjectile projectile(*incident);
289 hfs = bert->ApplyYourself(projectile, target);
290 } else {
291 // convert incident gamma to a pi0
293 G4double piKE = incident->GetTotalEnergy() - piMass;
294 G4double piMom = std::sqrt(piKE*(piKE + 2*piMass) );
295 G4ThreeVector piMomentum(incident->GetMomentumDirection() );
296 piMomentum *= piMom;
297 G4DynamicParticle theHadron(G4PionZero::PionZero(), piMomentum);
298 G4HadProjectile projectile(theHadron);
299 hfs = ftfp->ApplyYourself(projectile, target);
300 }
301
302 delete incident;
303
304 // Assign the creator model ID to the secondaries
305 for ( size_t i = 0; i < hfs->GetNumberOfSecondaries(); ++i ) {
307 }
308
309 // Copy secondaries from sub-model to model
311}
static constexpr double GeV
Definition: G4SIunits.hh:203
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
const G4ThreeVector & GetMomentumDirection() const
G4double GetTotalEnergy() const
void AddSecondaries(const std::vector< G4HadSecondary > &addSecs)
std::size_t GetNumberOfSecondaries() const
G4HadSecondary * GetSecondary(size_t i)
void SetCreatorModelID(G4int id)
static G4PionZero * PionZero()
Definition: G4PionZero.cc:107
G4HadFinalState * ApplyYourself(const G4HadProjectile &thePrimary, G4Nucleus &theNucleus) override

References G4HadFinalState::AddSecondaries(), G4CascadeInterface::ApplyYourself(), G4TheoFSGenerator::ApplyYourself(), bert, ftfp, G4DynamicParticle::GetMomentumDirection(), G4HadFinalState::GetNumberOfSecondaries(), G4ParticleDefinition::GetPDGMass(), G4HadFinalState::GetSecondary(), G4DynamicParticle::GetTotalEnergy(), GeV, G4PionZero::PionZero(), secID, G4HadSecondary::SetCreatorModelID(), and G4HadronicInteraction::theParticleChange.

Referenced by ApplyYourself().

◆ 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().

◆ 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(), 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

◆ GetVerboseLevel()

G4int G4HadronicInteraction::GetVerboseLevel ( ) const
inlineinherited

Definition at line 109 of file G4HadronicInteraction.hh.

References G4HadronicInteraction::verboseLevel.

◆ InitialiseModel()

void G4HadronicInteraction::InitialiseModel ( )
virtualinherited

◆ IsApplicable()

G4bool G4HadronicInteraction::IsApplicable ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtualinherited

◆ 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.

◆ MakeSamplingTable()

void G4MuonVDNuclearModel::MakeSamplingTable ( )
private

Definition at line 314 of file G4MuonVDNuclearModel.cc.

315{
316 G4int nbin;
318 G4double TotalEnergy;
319 G4double Maxep;
320 G4double CrossSection;
321
322 G4double c;
323 G4double y;
324 G4double ymin,ymax;
325 G4double dy,yy;
326 G4double dx,x;
327 G4double ep;
328
329 G4double AtomicNumber;
330 G4double AtomicWeight;
331
333
334 for (G4int iz = 0; iz < nzdat; ++iz) {
335 AtomicNumber = zdat[iz];
336 AtomicWeight = adat[iz]*(g/mole);
337
339 G4double pvv;
340
341 for (G4int it = 0; it < ntdat; ++it) {
342 KineticEnergy = tdat[it];
343 TotalEnergy = KineticEnergy + mumass;
344 Maxep = TotalEnergy - 0.5*proton_mass_c2;
345
346 CrossSection = 0.0;
347
348 // Calculate the differential cross section
349 // numerical integration in log .........
350 c = G4Log(Maxep/CutFixed);
351 ymin = -5.0;
352 ymax = 0.0;
353 dy = (ymax-ymin)/NBIN;
354
355 nbin=-1;
356
357 y = ymin - 0.5*dy;
358 yy = ymin - dy;
359 for (G4int i = 0; i < NBIN; ++i) {
360 y += dy;
361 x = G4Exp(y);
362 yy += dy;
363 dx = G4Exp(yy+dy)-G4Exp(yy);
364
365 ep = CutFixed*G4Exp(c*x);
366
367 CrossSection +=
368 ep*dx*muNucXS->ComputeDDMicroscopicCrossSection(KineticEnergy,
369 AtomicNumber,
370 AtomicWeight, ep);
371 if (nbin < NBIN) {
372 ++nbin;
373 pv->PutValue(nbin, it, CrossSection);
374 pv->PutX(nbin, y);
375 }
376 }
377 pv->PutX(NBIN, 0.);
378
379 if (CrossSection > 0.0) {
380 for (G4int ib = 0; ib <= nbin; ++ib) {
381 pvv = pv->GetValue(ib, it);
382 pvv = pvv/CrossSection;
383 pv->PutValue(ib, it, pvv);
384 }
385 }
386 } // loop on it
387
389 } // loop on iz
390
391 // G4cout << " Kokoulin XS = "
392 // << muNucXS->ComputeDDMicroscopicCrossSection(1*GeV, 20.0,
393 // 40.0*g/mole, 0.3*GeV)/millibarn
394 // << G4endl;
395}
static constexpr double mole
Definition: G4SIunits.hh:279
static constexpr double g
Definition: G4SIunits.hh:168
void InitialiseForElement(G4int Z, G4PhysicsVector *v)
G4double ComputeDDMicroscopicCrossSection(G4double incidentKE, G4double Z, G4double A, G4double epsilon)
static const G4double adat[nzdat]

References adat, G4KokoulinMuonNuclearXS::ComputeDDMicroscopicCrossSection(), CutFixed, fElementData, g, G4Exp(), G4Log(), G4ParticleDefinition::GetPDGMass(), G4ElementData::InitialiseForElement(), mole, muNucXS, G4MuonMinus::MuonMinus(), NBIN, ntdat, nzdat, source.hepunit::proton_mass_c2, tdat, and zdat.

Referenced by G4MuonVDNuclearModel().

◆ ModelDescription()

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

Reimplemented from G4HadronicInteraction.

Definition at line 397 of file G4MuonVDNuclearModel.cc.

398{
399 outFile << "G4MuonVDNuclearModel handles the inelastic scattering\n"
400 << "of mu- and mu+ from nuclei using the equivalent photon\n"
401 << "approximation in which the incoming lepton generates a\n"
402 << "virtual photon at the electromagnetic vertex, and the\n"
403 << "virtual photon is converted to a real photon. At low\n"
404 << "energies, the photon interacts directly with the nucleus\n"
405 << "using the Bertini cascade. At high energies the photon\n"
406 << "is converted to a pi0 which interacts using the FTFP\n"
407 << "model. The muon-nuclear cross sections of R. Kokoulin \n"
408 << "are used to generate the virtual photon spectrum\n";
409}

◆ operator!=()

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

◆ operator=()

G4MuonVDNuclearModel & G4MuonVDNuclearModel::operator= ( const G4MuonVDNuclearModel right)
privatedelete

◆ operator==()

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

◆ SampleInvariantT()

G4double G4HadronicInteraction::SampleInvariantT ( const G4ParticleDefinition p,
G4double  plab,
G4int  Z,
G4int  A 
)
virtualinherited

◆ 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().

◆ 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(), G4NeutrinoElectronCcModel::G4NeutrinoElectronCcModel(), G4NeutrinoElectronNcModel::G4NeutrinoElectronNcModel(), G4NeutrinoNucleusModel::G4NeutrinoNucleusModel(), G4NeutronElectronElModel::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(), G4NeutrinoElectronCcModel::G4NeutrinoElectronCcModel(), G4NeutrinoElectronNcModel::G4NeutrinoElectronNcModel(), G4NeutrinoNucleusModel::G4NeutrinoNucleusModel(), G4NeutronElectronElModel::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

Field Documentation

◆ adat

const G4double G4MuonVDNuclearModel::adat = {1.01,9.01,26.98,63.55,238.03}
staticprivate

Definition at line 79 of file G4MuonVDNuclearModel.hh.

Referenced by MakeSamplingTable().

◆ bert

G4CascadeInterface* G4MuonVDNuclearModel::bert
private

Definition at line 92 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateHadronicVertex(), and G4MuonVDNuclearModel().

◆ CutFixed

G4double G4MuonVDNuclearModel::CutFixed
private

◆ epCheckLevels

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

◆ fElementData

G4ElementData * G4MuonVDNuclearModel::fElementData = nullptr
staticprivate

◆ ftfp

G4TheoFSGenerator* G4MuonVDNuclearModel::ftfp
private

Definition at line 89 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateHadronicVertex(), and G4MuonVDNuclearModel().

◆ isBlocked

G4bool G4HadronicInteraction::isBlocked
protectedinherited

◆ isMaster

G4bool G4MuonVDNuclearModel::isMaster
private

Definition at line 85 of file G4MuonVDNuclearModel.hh.

Referenced by G4MuonVDNuclearModel(), and ~G4MuonVDNuclearModel().

◆ muNucXS

G4KokoulinMuonNuclearXS* G4MuonVDNuclearModel::muNucXS
private

Definition at line 87 of file G4MuonVDNuclearModel.hh.

Referenced by G4MuonVDNuclearModel(), and MakeSamplingTable().

◆ NBIN

const G4int G4MuonVDNuclearModel::NBIN = 800
staticprivate

Definition at line 74 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateEMVertex(), and MakeSamplingTable().

◆ ntdat

const G4int G4MuonVDNuclearModel::ntdat = 73
staticprivate

Definition at line 76 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateEMVertex(), and MakeSamplingTable().

◆ nzdat

const G4int G4MuonVDNuclearModel::nzdat = 5
staticprivate

Definition at line 75 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateEMVertex(), and MakeSamplingTable().

◆ recoilEnergyThreshold

G4double G4HadronicInteraction::recoilEnergyThreshold
privateinherited

◆ registry

G4HadronicInteractionRegistry* G4HadronicInteraction::registry
privateinherited

◆ secID

G4int G4MuonVDNuclearModel::secID
private

Definition at line 94 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateHadronicVertex(), and G4MuonVDNuclearModel().

◆ tdat

const G4double G4MuonVDNuclearModel::tdat
staticprivate
Initial value:
= {
1.e3,2.e3,3.e3,4.e3,5.e3,6.e3,7.e3,8.e3,9.e3,
1.e4,2.e4,3.e4,4.e4,5.e4,6.e4,7.e4,8.e4,9.e4,
1.e5,2.e5,3.e5,4.e5,5.e5,6.e5,7.e5,8.e5,9.e5,
1.e6,2.e6,3.e6,4.e6,5.e6,6.e6,7.e6,8.e6,9.e6,
1.e7,2.e7,3.e7,4.e7,5.e7,6.e7,7.e7,8.e7,9.e7,
1.e8,2.e8,3.e8,4.e8,5.e8,6.e8,7.e8,8.e8,9.e8,
1.e9,2.e9,3.e9,4.e9,5.e9,6.e9,7.e9,8.e9,9.e9,
1.e10,2.e10,3.e10,4.e10,5.e10,6.e10,7.e10,8.e10,9.e10,1.e11}

Definition at line 80 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateEMVertex(), and MakeSamplingTable().

◆ theBlockedList

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

◆ theBlockedListElements

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

◆ theFragmentation

G4LundStringFragmentation* G4MuonVDNuclearModel::theFragmentation
private

Definition at line 90 of file G4MuonVDNuclearModel.hh.

Referenced by G4MuonVDNuclearModel(), and ~G4MuonVDNuclearModel().

◆ 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

◆ 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(), G4NeutronElectronElModel::ApplyYourself(), G4LFission::ApplyYourself(), G4ANuElNucleusCcModel::ApplyYourself(), G4ANuElNucleusNcModel::ApplyYourself(), G4ANuMuNucleusCcModel::ApplyYourself(), G4ANuMuNucleusNcModel::ApplyYourself(), 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(), CalculateEMVertex(), G4ElectroVDNuclearModel::CalculateHadronicVertex(), 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().

◆ theStringDecay

G4ExcitedStringDecay* G4MuonVDNuclearModel::theStringDecay
private

Definition at line 91 of file G4MuonVDNuclearModel.hh.

Referenced by G4MuonVDNuclearModel(), and ~G4MuonVDNuclearModel().

◆ 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().

◆ zdat

const G4int G4MuonVDNuclearModel::zdat = {1, 4, 13, 29, 92}
staticprivate

Definition at line 78 of file G4MuonVDNuclearModel.hh.

Referenced by CalculateEMVertex(), and MakeSamplingTable().


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