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

INCL++ intra-nuclear cascade. More...

#include <G4INCLXXInterface.hh>

Inheritance diagram for G4INCLXXInterface:
G4VIntraNuclearTransportModel G4HadronicInteraction

Public Member Functions

void ActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Material *aMaterial)
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
void DeActivateFor (const G4Element *anElement)
 
void DeActivateFor (const G4Material *aMaterial)
 
void DeleteModel ()
 
 G4INCLXXInterface (G4VPreCompoundModel *const aPreCompound=0)
 
G4String const & GetDeExcitationModelName () const
 
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 G4VIntraNuclearTransportModel &right) const =delete
 
G4bool operator!= (G4INCLXXInterface &right)
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4bool operator== (const G4VIntraNuclearTransportModel &right) const =delete
 
G4bool operator== (G4INCLXXInterface &right)
 
G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
 
virtual void PropagateModelDescription (std::ostream &outFile) const
 
virtual G4ReactionProductVectorPropagateNuclNucl (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4V3DNucleus *theProjectileNucleus)
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
void Set3DNucleus (G4V3DNucleus *const value)
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
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 SetPrimaryProjectile (const G4HadProjectile &aPrimary)
 
void SetRecoilEnergyThreshold (G4double val)
 
void SetVerboseLevel (G4int value)
 
 ~G4INCLXXInterface ()
 

Protected Member Functions

void Block ()
 
G4V3DNucleusGet3DNucleus () const
 
G4VPreCompoundModelGetDeExcitation () const
 
const G4HadProjectileGetPrimaryProjectile () const
 
G4bool IsBlocked () const
 
void SetModelName (const G4String &nam)
 

Protected Attributes

G4bool isBlocked
 
G4V3DNucleusthe3DNucleus
 
G4VPreCompoundModeltheDeExcitation
 
G4double theMaxEnergy
 
G4double theMinEnergy
 
G4HadFinalState theParticleChange
 
const G4HadProjectilethePrimaryProjectile
 
G4String theTransportModelName
 
G4int verboseLevel
 

Private Member Functions

G4bool AccurateProjectile (const G4HadProjectile &aTrack, const G4Nucleus &theTargetNucleus) const
 
 G4INCLXXInterface (const G4INCLXXInterface &rhs)
 Dummy copy constructor to shut up Coverity warnings. More...
 
G4INCLXXInterfaceoperator= (G4INCLXXInterface const &rhs)
 Dummy assignment operator to shut up Coverity warnings. More...
 
G4double remnant4MomentumScaling (G4double mass, G4double kineticE, G4double px, G4double py, G4double pz) const
 Rescale remnant momentum if necessary. More...
 
G4DynamicParticletoG4Particle (G4int A, G4int Z, G4int S, G4int PDGCode, G4double kinE, G4double px, G4double py, G4double pz) const
 Convert an INCL particle to a G4DynamicParticle. More...
 
G4ParticleDefinitiontoG4ParticleDefinition (G4int A, G4int Z, G4int S, G4int PDGCode) const
 Convert A, Z and S to a G4ParticleDefinition. More...
 
G4double toINCLKineticEnergy (G4HadProjectile const &) const
 Convert G4HadProjectile to corresponding INCL particle kinetic energy. More...
 
G4INCL::ParticleSpecies toINCLParticleSpecies (G4HadProjectile const &) const
 Convert G4HadProjectile to corresponding INCL particle species. More...
 
G4INCL::ParticleType toINCLParticleType (G4ParticleDefinition const *const) const
 Convert G4ParticleDefinition to corresponding INCL particle type. More...
 

Private Attributes

G4bool complainedAboutBackupModel
 
G4bool complainedAboutPreCompound
 
G4bool dumpRemnantInfo
 
std::pair< G4double, G4doubleepCheckLevels
 
G4double recoilEnergyThreshold
 
G4HadronicInteractionRegistryregistry
 
G4int secID
 
G4HadronicInteractiontheBackupModel
 
G4HadronicInteractiontheBackupModelNucleon
 
std::vector< const G4Material * > theBlockedList
 
std::vector< const G4Element * > theBlockedListElements
 
G4INCL::INCLtheINCLModel
 
G4FissionProbabilitytheINCLXXFissionProbability
 
G4VLevelDensityParametertheINCLXXLevelDensity
 
G4INCLXXInterfaceStore *const theInterfaceStore
 
G4IonTable *const theIonTable
 
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
 
G4VPreCompoundModelthePreCompoundModel
 
G4HadFinalState theResult
 
G4INCLXXVInterfaceTallytheTally
 

Detailed Description

INCL++ intra-nuclear cascade.

Interface for INCL++. This interface handles basic hadron bullet particles (protons, neutrons, pions), as well as light ions.

Example usage in case of protons:

inclModel -> SetMinEnergy(0.0 * MeV); // Set the energy limits
inclModel -> SetMaxEnergy(3.0 * GeV);
G4HadronInelasticProcess* protonInelasticProcess = new G4HadronInelasticProcess( "protonInelastic", G4Proton::Definition() );
G4VCrossSectionDataSet* protonInelasticCrossSection = new G4BGGNucleonInelasticXS( G4Proton::Proton() );
protonInelasticProcess -> RegisterMe(inclModel);
protonInelasticProcess -> AddDataSet(protonInelasticCrossSection);
particle = G4Proton::Proton();
processManager = particle -> GetProcessManager();
processManager -> AddDiscreteProcess(protonInelasticProcess);
static constexpr double GeV
Definition: G4SIunits.hh:203
static constexpr double MeV
Definition: G4SIunits.hh:200
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
INCL++ intra-nuclear cascade.
G4INCLXXInterface(G4VPreCompoundModel *const aPreCompound=0)
static G4Proton * Definition()
Definition: G4Proton.cc:48
static G4Proton * Proton()
Definition: G4Proton.cc:92

The same setup procedure is needed for neutron, pion and generic-ion inelastic processes as well.

Definition at line 103 of file G4INCLXXInterface.hh.

Constructor & Destructor Documentation

◆ G4INCLXXInterface() [1/2]

G4INCLXXInterface::G4INCLXXInterface ( G4VPreCompoundModel *const  aPreCompound = 0)

Definition at line 68 of file G4INCLXXInterface.cc.

68 :
70 theINCLModel(NULL),
71 thePreCompoundModel(aPreCompound),
73 theTally(NULL),
79 secID(-1)
80{
86 }
87
88 // Use the environment variable G4INCLXX_NO_DE_EXCITATION to disable de-excitation
89 if(std::getenv("G4INCLXX_NO_DE_EXCITATION")) {
90 G4String message = "de-excitation is completely disabled!";
93 } else {
96 theDeExcitation = static_cast<G4VPreCompoundModel*>(p);
98
99 // set the fission parameters for G4ExcitationHandler
100 G4VEvaporationChannel * const theFissionChannel =
102 G4CompetitiveFission * const theFissionChannelCast = dynamic_cast<G4CompetitiveFission *>(theFissionChannel);
103 if(theFissionChannelCast) {
105 theFissionChannelCast->SetLevelDensityParameter(theINCLXXLevelDensity);
108 theFissionChannelCast->SetEmissionStrategy(theINCLXXFissionProbability);
109 theInterfaceStore->EmitBigWarning("INCL++/G4ExcitationHandler uses its own level-density parameter for fission");
110 } else {
111 theInterfaceStore->EmitBigWarning("INCL++/G4ExcitationHandler could not use its own level-density parameter for fission");
112 }
113 }
114
115 // use the envvar G4INCLXX_DUMP_REMNANT to dump information about the
116 // remnants on stdout
117 if(std::getenv("G4INCLXX_DUMP_REMNANT"))
118 dumpRemnantInfo = true;
119 else
120 dumpRemnantInfo = false;
121
124 secID = G4PhysicsModelCatalog::GetModelID( "model_INCLXXCascade" );
125}
void SetEmissionStrategy(G4VEmissionProbability *aFissionProb)
void SetLevelDensityParameter(G4VLevelDensityParameter *aLevelDensity)
G4VEvaporation * GetEvaporation()
Revised level-density parameter for fission after INCL++.
void SetFissionLevelDensityParameter(G4VLevelDensityParameter *aLevelDensity)
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
void EmitWarning(const G4String &message)
Emit a warning to G4cout.
static G4INCLXXInterfaceStore * GetInstance()
Get the singleton instance.
void EmitBigWarning(const G4String &message) const
Emit a BIG warning to G4cout.
G4INCLXXInterfaceStore *const theInterfaceStore
G4FissionProbability * theINCLXXFissionProbability
G4HadronicInteraction * theBackupModelNucleon
G4VLevelDensityParameter * theINCLXXLevelDensity
G4HadronicInteraction * theBackupModel
G4INCLXXVInterfaceTally * theTally
G4IonTable *const theIonTable
G4VPreCompoundModel * thePreCompoundModel
G4INCL::INCL * theINCLModel
static G4IonTable * GetIonTable()
Definition: G4IonTable.cc:170
static G4int GetModelID(const G4int modelIndex)
G4VEvaporationChannel * GetFissionChannel()
G4VIntraNuclearTransportModel(const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
G4ExcitationHandler * GetExcitationHandler() const

References dumpRemnantInfo, G4INCLXXInterfaceStore::EmitBigWarning(), G4INCLXXInterfaceStore::EmitWarning(), G4HadronicInteractionRegistry::FindModel(), G4ExcitationHandler::GetEvaporation(), G4VPreCompoundModel::GetExcitationHandler(), G4VEvaporation::GetFissionChannel(), G4PhysicsModelCatalog::GetModelID(), G4HadronicInteractionRegistry::Instance(), secID, G4CompetitiveFission::SetEmissionStrategy(), G4FissionProbability::SetFissionLevelDensityParameter(), G4CompetitiveFission::SetLevelDensityParameter(), theBackupModel, theBackupModelNucleon, G4VIntraNuclearTransportModel::theDeExcitation, theINCLXXFissionProbability, theINCLXXLevelDensity, theInterfaceStore, and thePreCompoundModel.

◆ ~G4INCLXXInterface()

G4INCLXXInterface::~G4INCLXXInterface ( )

Definition at line 127 of file G4INCLXXInterface.cc.

128{
131}

References theINCLXXFissionProbability, and theINCLXXLevelDensity.

◆ G4INCLXXInterface() [2/2]

G4INCLXXInterface::G4INCLXXInterface ( const G4INCLXXInterface rhs)
private

Dummy copy constructor to shut up Coverity warnings.

Member Function Documentation

◆ AccurateProjectile()

G4bool G4INCLXXInterface::AccurateProjectile ( const G4HadProjectile aTrack,
const G4Nucleus theTargetNucleus 
) const
private

Definition at line 133 of file G4INCLXXInterface.cc.

133 {
134 // Use direct kinematics if the projectile is a nucleon or a pion
135 const G4ParticleDefinition *projectileDef = aTrack.GetDefinition();
136 if(std::abs(projectileDef->GetBaryonNumber()) < 2) // Every non-composite particle. abs()-> in case of anti-nucleus projectile
137 return false;
138
139 // Here all projectiles should be nuclei
140 const G4int pA = projectileDef->GetAtomicMass();
141 if(pA<=0) {
142 std::stringstream ss;
143 ss << "the model does not know how to handle a collision between a "
144 << projectileDef->GetParticleName() << " projectile and a Z="
145 << theNucleus.GetZ_asInt() << ", A=" << theNucleus.GetA_asInt();
147 return true;
148 }
149
150 // If either nucleus is a LCP (A<=4), run the collision as light on heavy
151 const G4int tA = theNucleus.GetA_asInt();
152 if(tA<=4 || pA<=4) {
153 if(pA<tA)
154 return false;
155 else
156 return true;
157 }
158
159 // If one of the nuclei is heavier than theMaxProjMassINCL, run the collision
160 // as light on heavy.
161 // Note that here we are sure that either the projectile or the target is
162 // smaller than theMaxProjMass; otherwise theBackupModel would have been
163 // called.
164 const G4int theMaxProjMassINCL = theInterfaceStore->GetMaxProjMassINCL();
165 if(pA > theMaxProjMassINCL)
166 return true;
167 else if(tA > theMaxProjMassINCL)
168 return false;
169 else
170 // In all other cases, use the global setting
172}
int G4int
Definition: G4Types.hh:85
const G4ParticleDefinition * GetDefinition() const
G4int GetMaxProjMassINCL() const
Getter for theMaxProjMassINCL.
G4bool GetAccurateProjectile() const
Getter for accurateProjectile.
G4int GetAtomicMass() const
const G4String & GetParticleName() const

References G4INCLXXInterfaceStore::EmitBigWarning(), G4Nucleus::GetA_asInt(), G4INCLXXInterfaceStore::GetAccurateProjectile(), G4ParticleDefinition::GetAtomicMass(), G4ParticleDefinition::GetBaryonNumber(), G4HadProjectile::GetDefinition(), G4INCLXXInterfaceStore::GetMaxProjMassINCL(), G4ParticleDefinition::GetParticleName(), G4Nucleus::GetZ_asInt(), and theInterfaceStore.

Referenced by ApplyYourself().

◆ ActivateFor() [1/2]

void G4HadronicInteraction::ActivateFor ( const G4Element anElement)
inlineinherited

◆ ActivateFor() [2/2]

void G4HadronicInteraction::ActivateFor ( const G4Material aMaterial)
inlineinherited

◆ ApplyYourself()

G4HadFinalState * G4INCLXXInterface::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus theNucleus 
)
virtual

Main method to apply the INCL physics model.

Parameters
aTrackthe projectile particle
theNucleustarget nucleus
Returns
the output of the INCL physics model

Reimplemented from G4HadronicInteraction.

Definition at line 174 of file G4INCLXXInterface.cc.

175{
176 G4ParticleDefinition const * const trackDefinition = aTrack.GetDefinition();
177 const G4bool isIonTrack = trackDefinition->GetParticleType()==G4GenericIon::GenericIon()->GetParticleType();
178 const G4int trackA = trackDefinition->GetAtomicMass();
179 const G4int trackZ = (G4int) trackDefinition->GetPDGCharge();
180 const G4int trackL = trackDefinition->GetNumberOfLambdasInHypernucleus();
181 const G4int nucleusA = theNucleus.GetA_asInt();
182 const G4int nucleusZ = theNucleus.GetZ_asInt();
183
184 // For reactions induced by weird projectiles (e.g. He2), bail out
185 if((isIonTrack && ((trackZ<=0 && trackL==0) || trackA<=trackZ)) ||
186 (nucleusA>1 && (nucleusZ<=0 || nucleusA<=nucleusZ))) {
191 return &theResult;
192 }
193
194 // For reactions on nucleons, use the backup model (without complaining)
195 if(trackA<=1 && nucleusA<=1) {
196 return theBackupModelNucleon->ApplyYourself(aTrack, theNucleus);
197 }
198
199 // For systems heavier than theMaxProjMassINCL, use another model (typically
200 // BIC)
201 const G4int theMaxProjMassINCL = theInterfaceStore->GetMaxProjMassINCL();
202 if(trackA > theMaxProjMassINCL && nucleusA > theMaxProjMassINCL) {
205 std::stringstream ss;
206 ss << "INCL++ refuses to handle reactions between nuclei with A>"
207 << theMaxProjMassINCL
208 << ". A backup model ("
210 << ") will be used instead.";
212 }
213 return theBackupModel->ApplyYourself(aTrack, theNucleus);
214 }
215
216 // For energies lower than cascadeMinEnergyPerNucleon, use PreCompound
217 const G4double cascadeMinEnergyPerNucleon = theInterfaceStore->GetCascadeMinEnergyPerNucleon();
218 const G4double trackKinE = aTrack.GetKineticEnergy();
219 if((trackDefinition==G4Neutron::NeutronDefinition() || trackDefinition==G4Proton::ProtonDefinition())
220 && trackKinE < cascadeMinEnergyPerNucleon) {
223 std::stringstream ss;
224 ss << "INCL++ refuses to handle nucleon-induced reactions below "
225 << cascadeMinEnergyPerNucleon / MeV
226 << " MeV. A PreCoumpound model ("
228 << ") will be used instead.";
230 }
231 return thePreCompoundModel->ApplyYourself(aTrack, theNucleus);
232 }
233
234 // Calculate the total four-momentum in the entrance channel
235 const G4double theNucleusMass = theIonTable->GetIonMass(nucleusZ, nucleusA);
236 const G4double theTrackMass = trackDefinition->GetPDGMass();
237 const G4double theTrackEnergy = trackKinE + theTrackMass;
238 const G4double theTrackMomentumAbs2 = theTrackEnergy*theTrackEnergy - theTrackMass*theTrackMass;
239 const G4double theTrackMomentumAbs = ((theTrackMomentumAbs2>0.0) ? std::sqrt(theTrackMomentumAbs2) : 0.0);
240 const G4ThreeVector theTrackMomentum = aTrack.Get4Momentum().getV().unit() * theTrackMomentumAbs;
241
242 G4LorentzVector goodTrack4Momentum(theTrackMomentum, theTrackEnergy);
243 G4LorentzVector fourMomentumIn;
244 fourMomentumIn.setE(theTrackEnergy + theNucleusMass);
245 fourMomentumIn.setVect(theTrackMomentum);
246
247 // Check if inverse kinematics should be used
248 const G4bool inverseKinematics = AccurateProjectile(aTrack, theNucleus);
249
250 // If we are running in inverse kinematics, redefine aTrack and theNucleus
251 G4LorentzRotation *toInverseKinematics = NULL;
252 G4LorentzRotation *toDirectKinematics = NULL;
253 G4HadProjectile const *aProjectileTrack = &aTrack;
254 G4Nucleus *theTargetNucleus = &theNucleus;
255 if(inverseKinematics) {
256 G4ParticleDefinition *oldTargetDef = theIonTable->GetIon(nucleusZ, nucleusA, 0);
257 const G4ParticleDefinition *oldProjectileDef = trackDefinition;
258
259 if(oldProjectileDef != 0 && oldTargetDef != 0) {
260 const G4int newTargetA = oldProjectileDef->GetAtomicMass();
261 const G4int newTargetZ = oldProjectileDef->GetAtomicNumber();
262 const G4int newTargetL = oldProjectileDef->GetNumberOfLambdasInHypernucleus();
263 if(newTargetA > 0 && newTargetZ > 0) {
264 // This should give us the same energy per nucleon
265 theTargetNucleus = new G4Nucleus(newTargetA, newTargetZ, newTargetL);
266 toInverseKinematics = new G4LorentzRotation(goodTrack4Momentum.boostVector());
267 G4LorentzVector theProjectile4Momentum(0.0, 0.0, 0.0, theNucleusMass);
268 G4DynamicParticle swappedProjectileParticle(oldTargetDef, (*toInverseKinematics) * theProjectile4Momentum);
269 aProjectileTrack = new G4HadProjectile(swappedProjectileParticle);
270 } else {
271 G4String message = "badly defined target after swapping. Falling back to normal (non-swapped) mode.";
273 toInverseKinematics = new G4LorentzRotation;
274 }
275 } else {
276 G4String message = "oldProjectileDef or oldTargetDef was null";
278 toInverseKinematics = new G4LorentzRotation;
279 }
280 }
281
282 // INCL assumes the projectile particle is going in the direction of
283 // the Z-axis. Here we construct proper rotation to convert the
284 // momentum vectors of the outcoming particles to the original
285 // coordinate system.
286 G4LorentzVector projectileMomentum = aProjectileTrack->Get4Momentum();
287
288 // INCL++ assumes that the projectile is going in the direction of
289 // the z-axis. In principle, if the coordinate system used by G4
290 // hadronic framework is defined differently we need a rotation to
291 // transform the INCL++ reaction products to the appropriate
292 // frame. Please note that it isn't necessary to apply this
293 // transform to the projectile because when creating the INCL++
294 // projectile particle; \see{toINCLKineticEnergy} needs to use only the
295 // projectile energy (direction is simply assumed to be along z-axis).
297 toZ.rotateZ(-projectileMomentum.phi());
298 toZ.rotateY(-projectileMomentum.theta());
299 G4RotationMatrix toLabFrame3 = toZ.inverse();
300 G4LorentzRotation toLabFrame(toLabFrame3);
301 // However, it turns out that the projectile given to us by G4
302 // hadronic framework is already going in the direction of the
303 // z-axis so this rotation is actually unnecessary. Both toZ and
304 // toLabFrame turn out to be unit matrices as can be seen by
305 // uncommenting the folowing two lines:
306 // G4cout <<"toZ = " << toZ << G4endl;
307 // G4cout <<"toLabFrame = " << toLabFrame << G4endl;
308
309 theResult.Clear(); // Make sure the output data structure is clean.
311
312 std::list<G4Fragment> remnants;
313
314 const G4int maxTries = 200;
315 G4int nTries = 0;
316 // INCL can generate transparent events. However, this is meaningful
317 // only in the standalone code. In Geant4 we must "force" INCL to
318 // produce a valid cascade.
319 G4bool eventIsOK = false;
320 do {
321 const G4INCL::ParticleSpecies theSpecies = toINCLParticleSpecies(*aProjectileTrack);
322 const G4double kineticEnergy = toINCLKineticEnergy(*aProjectileTrack);
323
324 // The INCL model will be created at the first use
326
327 const G4INCL::EventInfo eventInfo = theINCLModel->processEvent(theSpecies, kineticEnergy,
328 theTargetNucleus->GetA_asInt(),
329 theTargetNucleus->GetZ_asInt(),
330 -theTargetNucleus->GetL()); // Strangeness has opposite sign
331 // eventIsOK = !eventInfo.transparent && nTries < maxTries; // of the number of Lambdas
332 eventIsOK = !eventInfo.transparent;
333 if(eventIsOK) {
334
335 // If the collision was run in inverse kinematics, prepare to boost back
336 // all the reaction products
337 if(inverseKinematics) {
338 toDirectKinematics = new G4LorentzRotation(toInverseKinematics->inverse());
339 }
340
341 G4LorentzVector fourMomentumOut;
342
343 for(G4int i = 0; i < eventInfo.nParticles; ++i) {
344 G4int A = eventInfo.A[i];
345 G4int Z = eventInfo.Z[i];
346 G4int S = eventInfo.S[i]; // Strangeness
347 G4int PDGCode = eventInfo.PDGCode[i];
348 // G4cout <<"INCL particle A = " << A << " Z = " << Z << " S = " << S << G4endl;
349 G4double kinE = eventInfo.EKin[i];
350 G4double px = eventInfo.px[i];
351 G4double py = eventInfo.py[i];
352 G4double pz = eventInfo.pz[i];
353 G4DynamicParticle *p = toG4Particle(A, Z, S, PDGCode, kinE, px, py, pz);
354 if(p != 0) {
355 G4LorentzVector momentum = p->Get4Momentum();
356
357 // Apply the toLabFrame rotation
358 momentum *= toLabFrame;
359 // Apply the inverse-kinematics boost
360 if(inverseKinematics) {
361 momentum *= *toDirectKinematics;
362 momentum.setVect(-momentum.vect());
363 }
364
365 // Set the four-momentum of the reaction products
366 p->Set4Momentum(momentum);
367 fourMomentumOut += momentum;
369
370 } else {
371 G4String message = "the model produced a particle that couldn't be converted to Geant4 particle.";
373 }
374 }
375
376 for(G4int i = 0; i < eventInfo.nRemnants; ++i) {
377 const G4int A = eventInfo.ARem[i];
378 const G4int Z = eventInfo.ZRem[i];
379 const G4int S = eventInfo.SRem[i];
380 // G4cout <<"INCL particle A = " << A << " Z = " << Z << " S= " << S << G4endl;
381 const G4double kinE = eventInfo.EKinRem[i];
382 const G4double px = eventInfo.pxRem[i];
383 const G4double py = eventInfo.pyRem[i];
384 const G4double pz = eventInfo.pzRem[i];
385 G4ThreeVector spin(
386 eventInfo.jxRem[i]*hbar_Planck,
387 eventInfo.jyRem[i]*hbar_Planck,
388 eventInfo.jzRem[i]*hbar_Planck
389 );
390 const G4double excitationE = eventInfo.EStarRem[i];
391 G4double nuclearMass = excitationE;
392 if ( S == 0 ) {
393 nuclearMass += G4NucleiProperties::GetNuclearMass(A, Z);
394 } else {
395 // Assumed that the opposite of the strangeness of the remnant gives the number of Lambdas inside it
396 nuclearMass += G4HyperNucleiProperties::GetNuclearMass(A, Z, std::abs(S));
397 }
398 const G4double scaling = remnant4MomentumScaling(nuclearMass, kinE, px, py, pz);
399 G4LorentzVector fourMomentum(scaling * px, scaling * py, scaling * pz,
400 nuclearMass + kinE);
401 if(std::abs(scaling - 1.0) > 0.01) {
402 std::stringstream ss;
403 ss << "momentum scaling = " << scaling
404 << "\n Lorentz vector = " << fourMomentum
405 << ")\n A = " << A << ", Z = " << Z << ", S = " << S
406 << "\n E* = " << excitationE << ", nuclearMass = " << nuclearMass
407 << "\n remnant i=" << i << ", nRemnants=" << eventInfo.nRemnants
408 << "\n Reaction was: " << aTrack.GetKineticEnergy()/MeV
409 << "-MeV " << trackDefinition->GetParticleName() << " + "
410 << theIonTable->GetIonName(theNucleus.GetZ_asInt(), theNucleus.GetA_asInt(), 0)
411 << ", in " << (inverseKinematics ? "inverse" : "direct") << " kinematics.";
413 }
414
415 // Apply the toLabFrame rotation
416 fourMomentum *= toLabFrame;
417 spin *= toLabFrame3;
418 // Apply the inverse-kinematics boost
419 if(inverseKinematics) {
420 fourMomentum *= *toDirectKinematics;
421 fourMomentum.setVect(-fourMomentum.vect());
422 }
423
424 fourMomentumOut += fourMomentum;
425 G4Fragment remnant(A, Z, std::abs(S), fourMomentum); // Assumed that -strangeness gives the number of Lambdas
426 remnant.SetAngularMomentum(spin);
427 remnant.SetCreatorModelID(secID);
428 if(dumpRemnantInfo) {
429 G4cerr << "G4INCLXX_DUMP_REMNANT: " << remnant << " spin: " << spin << G4endl;
430 }
431 remnants.push_back(remnant);
432 }
433
434 // Check four-momentum conservation
435 const G4LorentzVector violation4Momentum = fourMomentumOut - fourMomentumIn;
436 const G4double energyViolation = std::abs(violation4Momentum.e());
437 const G4double momentumViolation = violation4Momentum.rho();
438 if(energyViolation > G4INCLXXInterfaceStore::GetInstance()->GetConservationTolerance()) {
439 std::stringstream ss;
440 ss << "energy conservation violated by " << energyViolation/MeV << " MeV in "
441 << aTrack.GetKineticEnergy()/MeV << "-MeV " << trackDefinition->GetParticleName()
442 << " + " << theIonTable->GetIonName(theNucleus.GetZ_asInt(), theNucleus.GetA_asInt(), 0)
443 << " inelastic reaction, in " << (inverseKinematics ? "inverse" : "direct") << " kinematics. Will resample.";
445 eventIsOK = false;
446 const G4int nSecondaries = theResult.GetNumberOfSecondaries();
447 for(G4int j=0; j<nSecondaries; ++j)
451 remnants.clear();
452 } else if(momentumViolation > G4INCLXXInterfaceStore::GetInstance()->GetConservationTolerance()) {
453 std::stringstream ss;
454 ss << "momentum conservation violated by " << momentumViolation/MeV << " MeV in "
455 << aTrack.GetKineticEnergy()/MeV << "-MeV " << trackDefinition->GetParticleName()
456 << " + " << theIonTable->GetIonName(theNucleus.GetZ_asInt(), theNucleus.GetA_asInt(), 0)
457 << " inelastic reaction, in " << (inverseKinematics ? "inverse" : "direct") << " kinematics. Will resample.";
459 eventIsOK = false;
460 const G4int nSecondaries = theResult.GetNumberOfSecondaries();
461 for(G4int j=0; j<nSecondaries; ++j)
465 remnants.clear();
466 }
467 }
468 nTries++;
469 } while(!eventIsOK && nTries < maxTries); /* Loop checking, 10.07.2015, D.Mancusi */
470
471 // Clean up the objects that we created for the inverse kinematics
472 if(inverseKinematics) {
473 delete aProjectileTrack;
474 delete theTargetNucleus;
475 delete toInverseKinematics;
476 delete toDirectKinematics;
477 }
478
479 if(!eventIsOK) {
480 std::stringstream ss;
481 ss << "maximum number of tries exceeded for the proposed "
482 << aTrack.GetKineticEnergy()/MeV << "-MeV " << trackDefinition->GetParticleName()
483 << " + " << theIonTable->GetIonName(nucleusZ, nucleusA, 0)
484 << " inelastic reaction, in " << (inverseKinematics ? "inverse" : "direct") << " kinematics.";
489 return &theResult;
490 }
491
492 // De-excitation:
493
494 if(theDeExcitation != 0) {
495 for(std::list<G4Fragment>::iterator i = remnants.begin();
496 i != remnants.end(); ++i) {
497 G4ReactionProductVector *deExcitationResult = theDeExcitation->DeExcite((*i));
498
499 for(G4ReactionProductVector::iterator fragment = deExcitationResult->begin();
500 fragment != deExcitationResult->end(); ++fragment) {
501 const G4ParticleDefinition *def = (*fragment)->GetDefinition();
502 if(def != 0) {
503 G4DynamicParticle *theFragment = new G4DynamicParticle(def, (*fragment)->GetMomentum());
504 theResult.AddSecondary(theFragment, (*fragment)->GetCreatorModelID());
505 }
506 }
507
508 for(G4ReactionProductVector::iterator fragment = deExcitationResult->begin();
509 fragment != deExcitationResult->end(); ++fragment) {
510 delete (*fragment);
511 }
512 deExcitationResult->clear();
513 delete deExcitationResult;
514 }
515 }
516
517 remnants.clear();
518
520 theTally->Tally(aTrack, theNucleus, theResult);
521
522 return &theResult;
523}
G4double S(G4double temp)
@ isAlive
@ stopAndKill
CLHEP::HepLorentzRotation G4LorentzRotation
std::vector< G4ReactionProduct * > G4ReactionProductVector
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
const G4int Z[17]
const G4double A[17]
G4GLOB_DLL std::ostream G4cerr
#define G4endl
Definition: G4ios.hh:57
Hep3Vector unit() const
HepLorentzRotation inverse() const
double theta() const
Hep3Vector getV() const
Hep3Vector vect() const
void setVect(const Hep3Vector &)
HepRotation inverse() const
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:87
HepRotation & rotateY(double delta)
Definition: Rotation.cc:74
G4LorentzVector Get4Momentum() const
void Set4Momentum(const G4LorentzVector &momentum)
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:92
void SetStatusChange(G4HadFinalStateStatus aS)
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
std::size_t GetNumberOfSecondaries() const
void SetEnergyChange(G4double anEnergy)
G4HadSecondary * GetSecondary(size_t i)
void SetMomentumChange(const G4ThreeVector &aV)
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
G4DynamicParticle * GetParticle()
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
const G4String & GetModelName() const
static G4double GetNuclearMass(G4int A, G4int Z, G4int L)
G4INCL::INCL * GetINCLModel()
Get the cached INCL model engine.
G4double GetConservationTolerance() const
Getter for conservationTolerance.
G4double GetCascadeMinEnergyPerNucleon() const
Getter for cascadeMinEnergyPerNucleon.
G4INCLXXVInterfaceTally * GetTally() const
Getter for the interface tally.
G4HadFinalState theResult
G4double remnant4MomentumScaling(G4double mass, G4double kineticE, G4double px, G4double py, G4double pz) const
Rescale remnant momentum if necessary.
G4INCL::ParticleSpecies toINCLParticleSpecies(G4HadProjectile const &) const
Convert G4HadProjectile to corresponding INCL particle species.
G4double toINCLKineticEnergy(G4HadProjectile const &) const
Convert G4HadProjectile to corresponding INCL particle kinetic energy.
G4DynamicParticle * toG4Particle(G4int A, G4int Z, G4int S, G4int PDGCode, G4double kinE, G4double px, G4double py, G4double pz) const
Convert an INCL particle to a G4DynamicParticle.
G4bool AccurateProjectile(const G4HadProjectile &aTrack, const G4Nucleus &theTargetNucleus) const
virtual void Tally(G4HadProjectile const &aTrack, G4Nucleus const &theNucleus, G4HadFinalState const &result)=0
const EventInfo & processEvent()
const G4String & GetIonName(G4int Z, G4int A, G4int lvl=0) const
Definition: G4IonTable.cc:1229
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:522
G4double GetIonMass(G4int Z, G4int A, G4int nL=0, G4int lvl=0) const
Definition: G4IonTable.cc:1517
static G4Neutron * NeutronDefinition()
Definition: G4Neutron.cc:98
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4int GetA_asInt() const
Definition: G4Nucleus.hh:99
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:105
G4int GetL() const
Definition: G4Nucleus.hh:108
G4int GetAtomicNumber() const
const G4String & GetParticleType() const
G4double GetPDGCharge() const
G4int GetNumberOfLambdasInHypernucleus() const
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:87
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)=0
float hbar_Planck
Definition: hepunit.py:263
Short_t S[maxSizeParticles]
Particle strangeness number.
Short_t nRemnants
Number of remnants.
Bool_t transparent
True if the event is transparent.
Short_t A[maxSizeParticles]
Particle mass number.
Float_t EKinRem[maxSizeRemnants]
Remnant kinetic energy [MeV].
Float_t jzRem[maxSizeRemnants]
Remnant angular momentum, z component [ ].
Float_t EStarRem[maxSizeRemnants]
Remnant excitation energy [MeV].
Short_t Z[maxSizeParticles]
Particle charge number.
Float_t EKin[maxSizeParticles]
Particle kinetic energy [MeV].
Short_t nParticles
Number of particles in the final state.
Float_t jxRem[maxSizeRemnants]
Remnant angular momentum, x component [ ].
Float_t px[maxSizeParticles]
Particle momentum, x component [MeV/c].
Short_t SRem[maxSizeRemnants]
Remnant strangeness number.
Float_t pxRem[maxSizeRemnants]
Remnant momentum, x component [MeV/c].
Float_t pyRem[maxSizeRemnants]
Remnant momentum, y component [MeV/c].
Float_t pzRem[maxSizeRemnants]
Remnant momentum, z component [MeV/c].
Float_t pz[maxSizeParticles]
Particle momentum, z component [MeV/c].
Int_t PDGCode[maxSizeParticles]
PDG numbering of the particles.
Float_t jyRem[maxSizeRemnants]
Remnant angular momentum, y component [ ].
Float_t py[maxSizeParticles]
Particle momentum, y component [MeV/c].
Short_t ARem[maxSizeRemnants]
Remnant mass number.
Short_t ZRem[maxSizeRemnants]
Remnant charge number.

References A, G4INCL::EventInfo::A, AccurateProjectile(), G4HadFinalState::AddSecondary(), G4HadronicInteraction::ApplyYourself(), G4INCL::EventInfo::ARem, CLHEP::HepLorentzVector::boostVector(), G4HadFinalState::Clear(), complainedAboutBackupModel, complainedAboutPreCompound, G4VPreCompoundModel::DeExcite(), dumpRemnantInfo, CLHEP::HepLorentzVector::e(), G4INCL::EventInfo::EKin, G4INCL::EventInfo::EKinRem, G4INCLXXInterfaceStore::EmitBigWarning(), G4INCLXXInterfaceStore::EmitWarning(), G4INCL::EventInfo::EStarRem, G4cerr, G4endl, G4GenericIon::GenericIon(), G4DynamicParticle::Get4Momentum(), G4HadProjectile::Get4Momentum(), G4Nucleus::GetA_asInt(), G4ParticleDefinition::GetAtomicMass(), G4ParticleDefinition::GetAtomicNumber(), G4INCLXXInterfaceStore::GetCascadeMinEnergyPerNucleon(), G4INCLXXInterfaceStore::GetConservationTolerance(), G4HadProjectile::GetDefinition(), G4INCLXXInterfaceStore::GetINCLModel(), G4INCLXXInterfaceStore::GetInstance(), G4IonTable::GetIon(), G4IonTable::GetIonMass(), G4IonTable::GetIonName(), G4HadProjectile::GetKineticEnergy(), G4Nucleus::GetL(), G4INCLXXInterfaceStore::GetMaxProjMassINCL(), G4HadronicInteraction::GetModelName(), G4NucleiProperties::GetNuclearMass(), G4HyperNucleiProperties::GetNuclearMass(), G4ParticleDefinition::GetNumberOfLambdasInHypernucleus(), G4HadFinalState::GetNumberOfSecondaries(), G4HadSecondary::GetParticle(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetParticleType(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), G4HadFinalState::GetSecondary(), G4INCLXXInterfaceStore::GetTally(), CLHEP::HepLorentzVector::getV(), G4Nucleus::GetZ_asInt(), source.hepunit::hbar_Planck, CLHEP::HepLorentzRotation::inverse(), CLHEP::HepRotation::inverse(), isAlive, G4INCL::EventInfo::jxRem, G4INCL::EventInfo::jyRem, G4INCL::EventInfo::jzRem, MeV, G4Neutron::NeutronDefinition(), G4INCL::EventInfo::nParticles, G4INCL::EventInfo::nRemnants, G4INCL::EventInfo::PDGCode, CLHEP::HepLorentzVector::phi(), G4INCL::INCL::processEvent(), G4Proton::ProtonDefinition(), G4INCL::EventInfo::px, G4INCL::EventInfo::pxRem, G4INCL::EventInfo::py, G4INCL::EventInfo::pyRem, G4INCL::EventInfo::pz, G4INCL::EventInfo::pzRem, remnant4MomentumScaling(), CLHEP::HepLorentzVector::rho(), CLHEP::HepRotation::rotateY(), CLHEP::HepRotation::rotateZ(), S(), G4INCL::EventInfo::S, secID, G4DynamicParticle::Set4Momentum(), G4Fragment::SetAngularMomentum(), G4Fragment::SetCreatorModelID(), CLHEP::HepLorentzVector::setE(), G4HadFinalState::SetEnergyChange(), G4HadFinalState::SetMomentumChange(), G4HadFinalState::SetStatusChange(), CLHEP::HepLorentzVector::setVect(), G4INCL::EventInfo::SRem, stopAndKill, G4INCLXXVInterfaceTally::Tally(), theBackupModel, theBackupModelNucleon, G4VIntraNuclearTransportModel::theDeExcitation, theINCLModel, theInterfaceStore, theIonTable, thePreCompoundModel, theResult, CLHEP::HepLorentzVector::theta(), theTally, toG4Particle(), toINCLKineticEnergy(), toINCLParticleSpecies(), G4INCL::EventInfo::transparent, CLHEP::Hep3Vector::unit(), CLHEP::HepLorentzVector::vect(), Z, G4INCL::EventInfo::Z, and G4INCL::EventInfo::ZRem.

◆ Block()

void G4HadronicInteraction::Block ( )
inlineprotectedinherited

◆ BuildPhysicsTable()

void G4HadronicInteraction::BuildPhysicsTable ( const G4ParticleDefinition )
virtualinherited

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

◆ DeleteModel()

void G4INCLXXInterface::DeleteModel ( )
inline

Definition at line 128 of file G4INCLXXInterface.hh.

128 {
129 delete theINCLModel;
130 theINCLModel = NULL;
131 }

References theINCLModel.

◆ Get3DNucleus()

G4V3DNucleus * G4VIntraNuclearTransportModel::Get3DNucleus ( ) const
inlineprotectedinherited

◆ GetDeExcitation()

G4VPreCompoundModel * G4VIntraNuclearTransportModel::GetDeExcitation ( ) const
inlineprotectedinherited

◆ GetDeExcitationModelName()

G4String const & G4INCLXXInterface::GetDeExcitationModelName ( ) const

◆ 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 & G4VIntraNuclearTransportModel::GetModelName ( ) const
inlineinherited

◆ GetPrimaryProjectile()

const G4HadProjectile * G4VIntraNuclearTransportModel::GetPrimaryProjectile ( ) const
inlineprotectedinherited

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

◆ ModelDescription()

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

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 631 of file G4INCLXXInterface.cc.

631 {
632 outFile
633 << "The Li�ge Intranuclear Cascade (INCL++) is a model for reactions induced\n"
634 << "by nucleons, pions and light ion on any nucleus. The reaction is\n"
635 << "described as an avalanche of binary nucleon-nucleon collisions, which can\n"
636 << "lead to the emission of energetic particles and to the formation of an\n"
637 << "excited thermalised nucleus (remnant). The de-excitation of the remnant is\n"
638 << "outside the scope of INCL++ and is typically described by another model.\n\n"
639 << "INCL++ has been reasonably well tested for nucleon (~50 MeV to ~15 GeV),\n"
640 << "pion (idem) and light-ion projectiles (up to A=18, ~10A MeV to 1A GeV).\n"
641 << "Most tests involved target nuclei close to the stability valley, with\n"
642 << "numbers between 4 and 250.\n\n"
643 << "Reference: D. Mancusi et al., Phys. Rev. C90 (2014) 054602\n\n";
644}

◆ operator!=() [1/3]

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

◆ operator!=() [2/3]

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

◆ operator!=() [3/3]

G4bool G4INCLXXInterface::operator!= ( G4INCLXXInterface right)
inline

Definition at line 112 of file G4INCLXXInterface.hh.

112 {
113 return (this != &right);
114 }

◆ operator=()

G4INCLXXInterface & G4INCLXXInterface::operator= ( G4INCLXXInterface const &  rhs)
private

Dummy assignment operator to shut up Coverity warnings.

◆ operator==() [1/3]

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

◆ operator==() [2/3]

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

◆ operator==() [3/3]

G4bool G4INCLXXInterface::operator== ( G4INCLXXInterface right)
inline

Definition at line 108 of file G4INCLXXInterface.hh.

108 {
109 return (this == &right);
110 }

◆ Propagate()

G4ReactionProductVector * G4INCLXXInterface::Propagate ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus 
)
virtual

Implements G4VIntraNuclearTransportModel.

Definition at line 525 of file G4INCLXXInterface.cc.

525 {
526 return 0;
527}

◆ PropagateModelDescription()

void G4VIntraNuclearTransportModel::PropagateModelDescription ( std::ostream &  outFile) const
virtualinherited

Reimplemented in G4BinaryCascade, and G4GeneratorPrecompoundInterface.

Definition at line 58 of file G4VIntraNuclearTransportModel.cc.

59{
60 outFile << "G4VIntraNuclearTransportModel is abstract class, missing description.\n";
61}

Referenced by G4TheoFSGenerator::ModelDescription().

◆ PropagateNuclNucl()

G4ReactionProductVector * G4VIntraNuclearTransportModel::PropagateNuclNucl ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus,
G4V3DNucleus theProjectileNucleus 
)
virtualinherited

Reimplemented in G4GeneratorPrecompoundInterface.

Definition at line 64 of file G4VIntraNuclearTransportModel.cc.

66{
67 G4Exception("G4VIntraNuclearTransportModel::Propagate()","G4VINT02",
69 "Propagate method for nucleus-nucleus interactions not implemented");
70 return nullptr;
71}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35

References FatalException, and G4Exception().

Referenced by G4TheoFSGenerator::ApplyYourself().

◆ remnant4MomentumScaling()

G4double G4INCLXXInterface::remnant4MomentumScaling ( G4double  mass,
G4double  kineticE,
G4double  px,
G4double  py,
G4double  pz 
) const
private

Rescale remnant momentum if necessary.

Definition at line 618 of file G4INCLXXInterface.cc.

621 {
622 const G4double p2 = px*px + py*py + pz*pz;
623 if(p2 > 0.0) {
624 const G4double pnew2 = kineticE*kineticE + 2.0*kineticE*mass;
625 return std::sqrt(pnew2)/std::sqrt(p2);
626 } else {
627 return 1.0;
628 }
629}

Referenced by ApplyYourself().

◆ SampleInvariantT()

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

◆ Set3DNucleus()

void G4VIntraNuclearTransportModel::Set3DNucleus ( G4V3DNucleus *const  value)
inlineinherited

Definition at line 124 of file G4VIntraNuclearTransportModel.hh.

125{
126 delete the3DNucleus; the3DNucleus = value;
127}

References G4VIntraNuclearTransportModel::the3DNucleus.

◆ SetDeExcitation()

void G4VIntraNuclearTransportModel::SetDeExcitation ( G4VPreCompoundModel ptr)
inline

Definition at line 81 of file G4VIntraNuclearTransportModel.hh.

136{
137 // previous pre-compound model will be deleted at the end of job
138 theDeExcitation = value;
139}

Referenced by G4INCLXXInterfaceStore::UseAblaDeExcitation().

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

◆ SetPrimaryProjectile()

void G4VIntraNuclearTransportModel::SetPrimaryProjectile ( const G4HadProjectile aPrimary)
inlineinherited

Definition at line 148 of file G4VIntraNuclearTransportModel.hh.

149{
150 // NOTE: Previous pointer is NOT deleted: passed by reference, no ownership
151 thePrimaryProjectile = &aPrimary;
152}

References G4VIntraNuclearTransportModel::thePrimaryProjectile.

Referenced by G4TheoFSGenerator::ApplyYourself().

◆ SetRecoilEnergyThreshold()

void G4HadronicInteraction::SetRecoilEnergyThreshold ( G4double  val)
inlineinherited

◆ SetVerboseLevel()

void G4HadronicInteraction::SetVerboseLevel ( G4int  value)
inlineinherited

◆ toG4Particle()

G4DynamicParticle * G4INCLXXInterface::toG4Particle ( G4int  A,
G4int  Z,
G4int  S,
G4int  PDGCode,
G4double  kinE,
G4double  px,
G4double  py,
G4double  pz 
) const
private

Convert an INCL particle to a G4DynamicParticle.

Definition at line 604 of file G4INCLXXInterface.cc.

606 {
607 const G4ParticleDefinition *def = toG4ParticleDefinition(A, Z, S, PDGCode);
608 if(def == 0) { // Check if we have a valid particle definition
609 return 0;
610 }
611 const G4double energy = kinE * MeV;
612 const G4ThreeVector momentum(px, py, pz);
613 const G4ThreeVector momentumDirection = momentum.unit();
614 G4DynamicParticle *p = new G4DynamicParticle(def, momentumDirection, energy);
615 return p;
616}
G4ParticleDefinition * toG4ParticleDefinition(G4int A, G4int Z, G4int S, G4int PDGCode) const
Convert A, Z and S to a G4ParticleDefinition.
G4double energy(const ThreeVector &p, const G4double m)

References A, G4INCL::KinematicsUtils::energy(), MeV, S(), toG4ParticleDefinition(), CLHEP::Hep3Vector::unit(), and Z.

Referenced by ApplyYourself().

◆ toG4ParticleDefinition()

G4ParticleDefinition * G4INCLXXInterface::toG4ParticleDefinition ( G4int  A,
G4int  Z,
G4int  S,
G4int  PDGCode 
) const
private

Convert A, Z and S to a G4ParticleDefinition.

Definition at line 568 of file G4INCLXXInterface.cc.

568 {
569 if (PDGCode == 2212) { return G4Proton::Proton();
570 } else if(PDGCode == 2112) { return G4Neutron::Neutron();
571 } else if(PDGCode == 211) { return G4PionPlus::PionPlus();
572 } else if(PDGCode == 111) { return G4PionZero::PionZero();
573 } else if(PDGCode == -211) { return G4PionMinus::PionMinus();
574
575 } else if(PDGCode == 221) { return G4Eta::Eta();
576 } else if(PDGCode == 22) { return G4Gamma::Gamma();
577
578 } else if(PDGCode == 3122) { return G4Lambda::Lambda();
579 } else if(PDGCode == 3222) { return G4SigmaPlus::SigmaPlus();
580 } else if(PDGCode == 3212) { return G4SigmaZero::SigmaZero();
581 } else if(PDGCode == 3112) { return G4SigmaMinus::SigmaMinus();
582 } else if(PDGCode == 321) { return G4KaonPlus::KaonPlus();
583 } else if(PDGCode == -321) { return G4KaonMinus::KaonMinus();
584 } else if(PDGCode == 130) { return G4KaonZeroLong::KaonZeroLong();
585 } else if(PDGCode == 310) { return G4KaonZeroShort::KaonZeroShort();
586
587 } else if(PDGCode == 1002) { return G4Deuteron::Deuteron();
588 } else if(PDGCode == 1003) { return G4Triton::Triton();
589 } else if(PDGCode == 2003) { return G4He3::He3();
590 } else if(PDGCode == 2004) { return G4Alpha::Alpha();
591 } else if(S != 0) { // Assumed that -S gives the number of Lambdas
592 if (A == 3 && Z == 1 && S == -1 ) return G4HyperTriton::Definition();
593 if (A == 4 && Z == 1 && S == -1 ) return G4HyperH4::Definition();
594 if (A == 4 && Z == 2 && S == -1 ) return G4HyperAlpha::Definition();
595 if (A == 4 && Z == 1 && S == -2 ) return G4DoubleHyperH4::Definition();
596 if (A == 4 && Z == 0 && S == -2 ) return G4DoubleHyperDoubleNeutron::Definition();
597 if (A == 5 && Z == 2 && S == -1 ) return G4HyperHe5::Definition();
598 } else if(A > 0 && Z > 0 && A > Z) { // Returns ground state ion definition.
599 return theIonTable->GetIon(Z, A, 0);
600 }
601 return 0; // Error, unrecognized particle
602}
static G4Alpha * Alpha()
Definition: G4Alpha.cc:88
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:93
static G4DoubleHyperDoubleNeutron * Definition()
static G4DoubleHyperH4 * Definition()
static G4Eta * Eta()
Definition: G4Eta.cc:108
static G4Gamma * Gamma()
Definition: G4Gamma.cc:85
static G4He3 * He3()
Definition: G4He3.cc:93
static G4HyperAlpha * Definition()
Definition: G4HyperAlpha.cc:45
static G4HyperH4 * Definition()
Definition: G4HyperH4.cc:45
static G4HyperHe5 * Definition()
Definition: G4HyperHe5.cc:45
static G4HyperTriton * Definition()
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:112
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:112
static G4KaonZeroLong * KaonZeroLong()
static G4KaonZeroShort * KaonZeroShort()
static G4Lambda * Lambda()
Definition: G4Lambda.cc:107
static G4Neutron * Neutron()
Definition: G4Neutron.cc:103
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:97
static G4PionPlus * PionPlus()
Definition: G4PionPlus.cc:97
static G4PionZero * PionZero()
Definition: G4PionZero.cc:107
static G4SigmaMinus * SigmaMinus()
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:107
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:101
static G4Triton * Triton()
Definition: G4Triton.cc:93

References A, G4Alpha::Alpha(), G4DoubleHyperDoubleNeutron::Definition(), G4DoubleHyperH4::Definition(), G4HyperAlpha::Definition(), G4HyperH4::Definition(), G4HyperHe5::Definition(), G4HyperTriton::Definition(), G4Deuteron::Deuteron(), G4Eta::Eta(), G4Gamma::Gamma(), G4IonTable::GetIon(), G4He3::He3(), G4KaonMinus::KaonMinus(), G4KaonPlus::KaonPlus(), G4KaonZeroLong::KaonZeroLong(), G4KaonZeroShort::KaonZeroShort(), G4Lambda::Lambda(), G4Neutron::Neutron(), G4PionMinus::PionMinus(), G4PionPlus::PionPlus(), G4PionZero::PionZero(), G4Proton::Proton(), S(), G4SigmaMinus::SigmaMinus(), G4SigmaPlus::SigmaPlus(), G4SigmaZero::SigmaZero(), theIonTable, G4Triton::Triton(), and Z.

Referenced by toG4Particle().

◆ toINCLKineticEnergy()

G4double G4INCLXXInterface::toINCLKineticEnergy ( G4HadProjectile const &  aTrack) const
private

Convert G4HadProjectile to corresponding INCL particle kinetic energy.

Definition at line 564 of file G4INCLXXInterface.cc.

564 {
565 return aTrack.GetKineticEnergy();
566}

References G4HadProjectile::GetKineticEnergy().

Referenced by ApplyYourself().

◆ toINCLParticleSpecies()

G4INCL::ParticleSpecies G4INCLXXInterface::toINCLParticleSpecies ( G4HadProjectile const &  aTrack) const
private

Convert G4HadProjectile to corresponding INCL particle species.

Definition at line 550 of file G4INCLXXInterface.cc.

550 {
551 const G4ParticleDefinition *pdef = aTrack.GetDefinition();
552 const G4INCL::ParticleType theType = toINCLParticleType(pdef);
553 if(theType!=G4INCL::Composite)
554 return G4INCL::ParticleSpecies(theType);
555 else {
556 G4INCL::ParticleSpecies theSpecies;
557 theSpecies.theType=theType;
558 theSpecies.theA=pdef->GetAtomicMass();
559 theSpecies.theZ=pdef->GetAtomicNumber();
560 return theSpecies;
561 }
562}
G4INCL::ParticleType toINCLParticleType(G4ParticleDefinition const *const) const
Convert G4ParticleDefinition to corresponding INCL particle type.

References G4INCL::Composite, G4ParticleDefinition::GetAtomicMass(), G4ParticleDefinition::GetAtomicNumber(), G4HadProjectile::GetDefinition(), G4INCL::ParticleSpecies::theA, G4INCL::ParticleSpecies::theType, G4INCL::ParticleSpecies::theZ, and toINCLParticleType().

Referenced by ApplyYourself().

◆ toINCLParticleType()

G4INCL::ParticleType G4INCLXXInterface::toINCLParticleType ( G4ParticleDefinition const * const  pdef) const
private

Convert G4ParticleDefinition to corresponding INCL particle type.

Definition at line 529 of file G4INCLXXInterface.cc.

529 {
530 if( pdef == G4Proton::Proton()) return G4INCL::Proton;
531 else if(pdef == G4Neutron::Neutron()) return G4INCL::Neutron;
532 else if(pdef == G4PionPlus::PionPlus()) return G4INCL::PiPlus;
533 else if(pdef == G4PionMinus::PionMinus()) return G4INCL::PiMinus;
534 else if(pdef == G4PionZero::PionZero()) return G4INCL::PiZero;
535 else if(pdef == G4KaonPlus::KaonPlus()) return G4INCL::KPlus;
536 else if(pdef == G4KaonZero::KaonZero()) return G4INCL::KZero;
537 else if(pdef == G4KaonMinus::KaonMinus()) return G4INCL::KMinus;
538 else if(pdef == G4AntiKaonZero::AntiKaonZero()) return G4INCL::KZeroBar;
539 // For K0L & K0S we do not take into account K0/K0B oscillations
540 else if(pdef == G4KaonZeroLong::KaonZeroLong()) return G4UniformRand() < 0.5 ? G4INCL::KZeroBar : G4INCL::KZero;
541 else if(pdef == G4KaonZeroShort::KaonZeroShort()) return G4UniformRand() < 0.5 ? G4INCL::KZeroBar : G4INCL::KZero;
542 else if(pdef == G4Deuteron::Deuteron()) return G4INCL::Composite;
543 else if(pdef == G4Triton::Triton()) return G4INCL::Composite;
544 else if(pdef == G4He3::He3()) return G4INCL::Composite;
545 else if(pdef == G4Alpha::Alpha()) return G4INCL::Composite;
546 else if(pdef->GetParticleType() == G4GenericIon::GenericIon()->GetParticleType()) return G4INCL::Composite;
547 else return G4INCL::UnknownParticle;
548}
#define G4UniformRand()
Definition: Randomize.hh:52
static G4AntiKaonZero * AntiKaonZero()
static G4KaonZero * KaonZero()
Definition: G4KaonZero.cc:103

References G4Alpha::Alpha(), G4AntiKaonZero::AntiKaonZero(), G4INCL::Composite, G4Deuteron::Deuteron(), G4UniformRand, G4GenericIon::GenericIon(), G4ParticleDefinition::GetParticleType(), G4He3::He3(), G4KaonMinus::KaonMinus(), G4KaonPlus::KaonPlus(), G4KaonZero::KaonZero(), G4KaonZeroLong::KaonZeroLong(), G4KaonZeroShort::KaonZeroShort(), G4INCL::KMinus, G4INCL::KPlus, G4INCL::KZero, G4INCL::KZeroBar, G4Neutron::Neutron(), G4INCL::Neutron, G4INCL::PiMinus, G4PionMinus::PionMinus(), G4PionPlus::PionPlus(), G4PionZero::PionZero(), G4INCL::PiPlus, G4INCL::PiZero, G4Proton::Proton(), G4INCL::Proton, G4Triton::Triton(), and G4INCL::UnknownParticle.

Referenced by toINCLParticleSpecies().

Field Documentation

◆ complainedAboutBackupModel

G4bool G4INCLXXInterface::complainedAboutBackupModel
private

Definition at line 178 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself().

◆ complainedAboutPreCompound

G4bool G4INCLXXInterface::complainedAboutPreCompound
private

Definition at line 179 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself().

◆ dumpRemnantInfo

G4bool G4INCLXXInterface::dumpRemnantInfo
private

Definition at line 183 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and G4INCLXXInterface().

◆ epCheckLevels

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

◆ isBlocked

G4bool G4HadronicInteraction::isBlocked
protectedinherited

◆ recoilEnergyThreshold

G4double G4HadronicInteraction::recoilEnergyThreshold
privateinherited

◆ registry

G4HadronicInteractionRegistry* G4HadronicInteraction::registry
privateinherited

◆ secID

G4int G4INCLXXInterface::secID
private

Definition at line 188 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and G4INCLXXInterface().

◆ the3DNucleus

G4V3DNucleus* G4VIntraNuclearTransportModel::the3DNucleus
protectedinherited

◆ theBackupModel

G4HadronicInteraction* G4INCLXXInterface::theBackupModel
private

Definition at line 172 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and G4INCLXXInterface().

◆ theBackupModelNucleon

G4HadronicInteraction* G4INCLXXInterface::theBackupModelNucleon
private

Definition at line 173 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and G4INCLXXInterface().

◆ theBlockedList

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

◆ theBlockedListElements

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

◆ theDeExcitation

G4VPreCompoundModel* G4VIntraNuclearTransportModel::theDeExcitation
protectedinherited

◆ theINCLModel

G4INCL::INCL* G4INCLXXInterface::theINCLModel
private

Definition at line 166 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and DeleteModel().

◆ theINCLXXFissionProbability

G4FissionProbability* G4INCLXXInterface::theINCLXXFissionProbability
private

Definition at line 186 of file G4INCLXXInterface.hh.

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

◆ theINCLXXLevelDensity

G4VLevelDensityParameter* G4INCLXXInterface::theINCLXXLevelDensity
private

Definition at line 185 of file G4INCLXXInterface.hh.

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

◆ theInterfaceStore

G4INCLXXInterfaceStore* const G4INCLXXInterface::theInterfaceStore
private

Definition at line 175 of file G4INCLXXInterface.hh.

Referenced by AccurateProjectile(), ApplyYourself(), and G4INCLXXInterface().

◆ theIonTable

G4IonTable* const G4INCLXXInterface::theIonTable
private

Definition at line 181 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and toG4ParticleDefinition().

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

◆ thePreCompoundModel

G4VPreCompoundModel* G4INCLXXInterface::thePreCompoundModel
private

Definition at line 168 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself(), and G4INCLXXInterface().

◆ thePrimaryProjectile

const G4HadProjectile* G4VIntraNuclearTransportModel::thePrimaryProjectile
protectedinherited

◆ theResult

G4HadFinalState G4INCLXXInterface::theResult
private

Definition at line 170 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself().

◆ theTally

G4INCLXXVInterfaceTally* G4INCLXXInterface::theTally
private

Definition at line 176 of file G4INCLXXInterface.hh.

Referenced by ApplyYourself().

◆ theTransportModelName

G4String G4VIntraNuclearTransportModel::theTransportModelName
protectedinherited

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