Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4INCLXXInterface Class Reference

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

#include <G4INCLXXInterface.hh>

+ Inheritance diagram for G4INCLXXInterface:

Public Member Functions

 G4INCLXXInterface (G4VPreCompoundModel *const aPreCompound=0)
 
 ~G4INCLXXInterface ()
 
G4bool operator== (G4INCLXXInterface &right)
 
G4bool operator!= (G4INCLXXInterface &right)
 
G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
void DeleteModel ()
 
virtual void ModelDescription (std::ostream &outFile) const
 
G4String const & GetDeExcitationModelName () const
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
- Public Member Functions inherited from G4VIntraNuclearTransportModel
 G4VIntraNuclearTransportModel (const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
 
virtual ~G4VIntraNuclearTransportModel ()
 
virtual G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)=0
 
virtual G4ReactionProductVectorPropagateNuclNucl (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4V3DNucleus *theProjectileNucleus)
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
void Set3DNucleus (G4V3DNucleus *const value)
 
void SetPrimaryProjectile (const G4HadProjectile &aPrimary)
 
const G4StringGetModelName () const
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void PropagateModelDescription (std::ostream &outFile) const
 
 G4VIntraNuclearTransportModel (const G4VIntraNuclearTransportModel &right)=delete
 
const G4VIntraNuclearTransportModeloperator= (const G4VIntraNuclearTransportModel &right)=delete
 
G4bool operator== (const G4VIntraNuclearTransportModel &right) const =delete
 
G4bool operator!= (const G4VIntraNuclearTransportModel &right) const =delete
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair< G4double, G4doubleGetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 
 G4HadronicInteraction (const G4HadronicInteraction &right)=delete
 
const G4HadronicInteractionoperator= (const G4HadronicInteraction &right)=delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 

Additional Inherited Members

- Protected Member Functions inherited from G4VIntraNuclearTransportModel
G4V3DNucleusGet3DNucleus () const
 
G4VPreCompoundModelGetDeExcitation () const
 
const G4HadProjectileGetPrimaryProjectile () const
 
- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4VIntraNuclearTransportModel
G4String theTransportModelName
 
G4V3DNucleusthe3DNucleus
 
G4VPreCompoundModeltheDeExcitation
 
const G4HadProjectilethePrimaryProjectile
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

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);
G4ProtonInelasticProcess* protonInelasticProcess = new G4ProtonInelasticProcess();
G4ProtonInelasticCrossSection* protonInelasticCrossSection = new G4ProtonInelasticCrossSection();
protonInelasticProcess -> RegisterMe(inclModel);
protonInelasticProcess -> AddDataSet(protonInelasticCrossSection);
particle = G4Proton::Proton();
processManager = particle -> GetProcessManager();
processManager -> AddDiscreteProcess(protonInelasticProcess);
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
INCL++ intra-nuclear cascade.
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()

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

Definition at line 59 of file G4INCLXXInterface.cc.

59 :
61 theINCLModel(NULL),
62 thePreCompoundModel(aPreCompound),
63 theInterfaceStore(G4INCLXXInterfaceStore::GetInstance()),
64 theTally(NULL),
65 complainedAboutBackupModel(false),
66 complainedAboutPreCompound(false),
67 theIonTable(G4IonTable::GetIonTable()),
68 theINCLXXLevelDensity(NULL),
69 theINCLXXFissionProbability(NULL)
70{
71 if(!thePreCompoundModel) {
74 thePreCompoundModel = static_cast<G4VPreCompoundModel*>(p);
75 if(!thePreCompoundModel) { thePreCompoundModel = new G4PreCompoundModel; }
76 }
77
78 // Use the environment variable G4INCLXX_NO_DE_EXCITATION to disable de-excitation
79 if(std::getenv("G4INCLXX_NO_DE_EXCITATION")) {
80 G4String message = "de-excitation is completely disabled!";
81 theInterfaceStore->EmitWarning(message);
83 } else {
86 theDeExcitation = static_cast<G4VPreCompoundModel*>(p);
88
89 // set the fission parameters for G4ExcitationHandler
90 G4VEvaporationChannel * const theFissionChannel =
92 G4CompetitiveFission * const theFissionChannelCast = dynamic_cast<G4CompetitiveFission *>(theFissionChannel);
93 if(theFissionChannelCast) {
94 theINCLXXLevelDensity = new G4FissionLevelDensityParameterINCLXX;
95 theFissionChannelCast->SetLevelDensityParameter(theINCLXXLevelDensity);
96 theINCLXXFissionProbability = new G4FissionProbability;
97 theINCLXXFissionProbability->SetFissionLevelDensityParameter(theINCLXXLevelDensity);
98 theFissionChannelCast->SetEmissionStrategy(theINCLXXFissionProbability);
99 theInterfaceStore->EmitBigWarning("INCL++/G4ExcitationHandler uses its own level-density parameter for fission");
100 } else {
101 theInterfaceStore->EmitBigWarning("INCL++/G4ExcitationHandler could not use its own level-density parameter for fission");
102 }
103 }
104
105 // use the envvar G4INCLXX_DUMP_REMNANT to dump information about the
106 // remnants on stdout
107 if(std::getenv("G4INCLXX_DUMP_REMNANT"))
108 dumpRemnantInfo = true;
109 else
110 dumpRemnantInfo = false;
111
112 theBackupModel = new G4BinaryLightIonReaction;
113 theBackupModelNucleon = new G4BinaryCascade;
114}
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.
static G4IonTable * GetIonTable()
Definition: G4IonTable.cc:170
G4VEvaporationChannel * GetFissionChannel()
G4ExcitationHandler * GetExcitationHandler() const

◆ ~G4INCLXXInterface()

G4INCLXXInterface::~G4INCLXXInterface ( )

Definition at line 116 of file G4INCLXXInterface.cc.

117{
118 delete theINCLXXLevelDensity;
119 delete theINCLXXFissionProbability;
120}

Member Function Documentation

◆ 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 163 of file G4INCLXXInterface.cc.

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

◆ DeleteModel()

void G4INCLXXInterface::DeleteModel ( )
inline

Definition at line 128 of file G4INCLXXInterface.hh.

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

◆ GetDeExcitationModelName()

G4String const & G4INCLXXInterface::GetDeExcitationModelName ( ) const

Definition at line 613 of file G4INCLXXInterface.cc.

613 {
615}

◆ ModelDescription()

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

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 598 of file G4INCLXXInterface.cc.

598 {
599 outFile
600 << "The Li�ge Intranuclear Cascade (INCL++) is a model for reactions induced\n"
601 << "by nucleons, pions and light ion on any nucleus. The reaction is\n"
602 << "described as an avalanche of binary nucleon-nucleon collisions, which can\n"
603 << "lead to the emission of energetic particles and to the formation of an\n"
604 << "excited thermalised nucleus (remnant). The de-excitation of the remnant is\n"
605 << "outside the scope of INCL++ and is typically described by another model.\n\n"
606 << "INCL++ has been reasonably well tested for nucleon (~50 MeV to ~15 GeV),\n"
607 << "pion (idem) and light-ion projectiles (up to A=18, ~10A MeV to 1A GeV).\n"
608 << "Most tests involved target nuclei close to the stability valley, with\n"
609 << "numbers between 4 and 250.\n\n"
610 << "Reference: D. Mancusi et al., Phys. Rev. C90 (2014) 054602\n\n";
611}

◆ operator!=()

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

Definition at line 112 of file G4INCLXXInterface.hh.

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

◆ operator==()

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 502 of file G4INCLXXInterface.cc.

502 {
503 return 0;
504}

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


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