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

#include <G4GeneratorPrecompoundInterface.hh>

Inheritance diagram for G4GeneratorPrecompoundInterface:
G4VIntraNuclearTransportModel G4HadronicInteraction

Public Member Functions

void ActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Material *aMaterial)
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
void DeActivateFor (const G4Element *anElement)
 
void DeActivateFor (const G4Material *aMaterial)
 
 G4GeneratorPrecompoundInterface (G4VPreCompoundModel *p=0)
 
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
 
void MakeCoalescence (G4KineticTrackVector *theSecondaries)
 
virtual void ModelDescription (std::ostream &outFile) const
 
G4bool operator!= (const G4HadronicInteraction &right) const =delete
 
G4bool operator!= (const G4VIntraNuclearTransportModel &right) const =delete
 
G4bool operator== (const G4HadronicInteraction &right) const =delete
 
G4bool operator== (const G4VIntraNuclearTransportModel &right) const =delete
 
virtual G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
 
virtual void PropagateModelDescription (std::ostream &) 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 SetCaptureThreshold (G4double)
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
void SetDeltaM (G4double)
 
void SetDeltaR (G4double)
 
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)
 
virtual ~G4GeneratorPrecompoundInterface ()
 

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

 G4GeneratorPrecompoundInterface (const G4GeneratorPrecompoundInterface &right)
 
G4bool operator!= (G4GeneratorPrecompoundInterface &right)
 
const G4GeneratorPrecompoundInterfaceoperator= (const G4GeneratorPrecompoundInterface &right)
 
G4bool operator== (G4GeneratorPrecompoundInterface &right)
 

Private Attributes

const G4ParticleDefinitionANTIdeuteron
 
const G4ParticleDefinitionANTIHe3
 
const G4ParticleDefinitionANTIHe4
 
const G4ParticleDefinitionANTIneutron
 
const G4ParticleDefinitionANTIproton
 
const G4ParticleDefinitionANTItriton
 
G4double CaptureThreshold
 
G4double DeltaM
 
G4double DeltaR
 
const G4ParticleDefinitiondeuteron
 
std::pair< G4double, G4doubleepCheckLevels
 
const G4ParticleDefinitionHe3
 
const G4ParticleDefinitionHe4
 
const G4ParticleDefinitionneutron
 
const G4ParticleDefinitionproton
 
G4double recoilEnergyThreshold
 
G4HadronicInteractionRegistryregistry
 
G4int secID
 
std::vector< const G4Material * > theBlockedList
 
std::vector< const G4Element * > theBlockedListElements
 
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
 
const G4ParticleDefinitiontriton
 

Detailed Description

Definition at line 59 of file G4GeneratorPrecompoundInterface.hh.

Constructor & Destructor Documentation

◆ G4GeneratorPrecompoundInterface() [1/2]

G4GeneratorPrecompoundInterface::G4GeneratorPrecompoundInterface ( G4VPreCompoundModel p = 0)

Definition at line 98 of file G4GeneratorPrecompoundInterface.cc.

99 : CaptureThreshold(70*MeV), DeltaM(5.0*MeV), DeltaR(0.0), secID(-1)
100{
103
106 He3 =G4He3::He3();
108
111
116
117 if(preModel) { SetDeExcitation(preModel); }
118 else {
121 G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(hadi);
122 if(!pre) { pre = new G4PreCompoundModel(); }
123 SetDeExcitation(pre);
124 }
125
127}
static constexpr double MeV
Definition: G4SIunits.hh:200
static G4Alpha * Alpha()
Definition: G4Alpha.cc:88
static G4AntiAlpha * AntiAlpha()
Definition: G4AntiAlpha.cc:88
static G4AntiDeuteron * AntiDeuteron()
static G4AntiHe3 * AntiHe3()
Definition: G4AntiHe3.cc:93
static G4AntiNeutron * AntiNeutron()
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:92
static G4AntiTriton * AntiTriton()
Definition: G4AntiTriton.cc:93
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:93
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
static G4He3 * He3()
Definition: G4He3.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:103
static G4int GetModelID(const G4int modelIndex)
static G4Proton * Proton()
Definition: G4Proton.cc:92
static G4Triton * Triton()
Definition: G4Triton.cc:93
void SetDeExcitation(G4VPreCompoundModel *ptr)

References G4Alpha::Alpha(), G4AntiAlpha::AntiAlpha(), G4AntiDeuteron::AntiDeuteron(), ANTIdeuteron, G4AntiHe3::AntiHe3(), ANTIHe3, ANTIHe4, G4AntiNeutron::AntiNeutron(), ANTIneutron, G4AntiProton::AntiProton(), ANTIproton, G4AntiTriton::AntiTriton(), ANTItriton, G4Deuteron::Deuteron(), deuteron, G4HadronicInteractionRegistry::FindModel(), G4PhysicsModelCatalog::GetModelID(), G4He3::He3(), He3, He4, G4HadronicInteractionRegistry::Instance(), G4Neutron::Neutron(), neutron, G4Proton::Proton(), proton, secID, G4VIntraNuclearTransportModel::SetDeExcitation(), G4Triton::Triton(), and triton.

◆ ~G4GeneratorPrecompoundInterface()

G4GeneratorPrecompoundInterface::~G4GeneratorPrecompoundInterface ( )
virtual

Definition at line 129 of file G4GeneratorPrecompoundInterface.cc.

130{
131}

◆ G4GeneratorPrecompoundInterface() [2/2]

G4GeneratorPrecompoundInterface::G4GeneratorPrecompoundInterface ( const G4GeneratorPrecompoundInterface right)
private

Member Function Documentation

◆ ActivateFor() [1/2]

void G4HadronicInteraction::ActivateFor ( const G4Element anElement)
inlineinherited

◆ ActivateFor() [2/2]

void G4HadronicInteraction::ActivateFor ( const G4Material aMaterial)
inlineinherited

◆ ApplyYourself()

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

Reimplemented from G4HadronicInteraction.

Definition at line 349 of file G4GeneratorPrecompoundInterface.cc.

351{
352 G4cout << "G4GeneratorPrecompoundInterface: ApplyYourself interface called stand-allone."
353 << G4endl;
354 G4cout << "This class is only a mediator between generator and precompound"<<G4endl;
355 G4cout << "Please remove from your physics list."<<G4endl;
356 throw G4HadronicException(__FILE__, __LINE__, "SEVERE: G4GeneratorPrecompoundInterface model interface called stand-allone.");
357 return new G4HadFinalState;
358}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout

References G4cout, and G4endl.

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

◆ Get3DNucleus()

G4V3DNucleus * G4VIntraNuclearTransportModel::Get3DNucleus ( ) const
inlineprotectedinherited

◆ GetDeExcitation()

G4VPreCompoundModel * G4VIntraNuclearTransportModel::GetDeExcitation ( ) const
inlineprotectedinherited

◆ 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
static constexpr double GeV
Definition: G4SIunits.hh:203

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.

◆ MakeCoalescence()

void G4GeneratorPrecompoundInterface::MakeCoalescence ( G4KineticTrackVector theSecondaries)

Definition at line 786 of file G4GeneratorPrecompoundInterface.cc.

786 {
787 // This method replaces pairs of proton-neutron - in the case of nuclei - or
788 // antiproton-antineutron - in the case of anti-nuclei - which are close in
789 // momentum, with, respectively, deuterons and anti-deuterons.
790 // Note that in the case of hypernuclei or anti-hypernuclei, lambdas or anti-lambdas
791 // are not considered for coalescence because hyper-deuteron or anti-hyper-deuteron
792 // are assumed not to exist.
793
794 if (!tracks) return;
795
796 G4double MassCut = deuteron->GetPDGMass() + DeltaM; // In MeV
797
798 for ( size_t i = 0; i < tracks->size(); ++i ) { // search for protons
799
800 G4KineticTrack* trackP = (*tracks)[i];
801 if ( ! trackP ) continue;
802 if (trackP->GetDefinition() != proton) continue;
803
804 G4LorentzVector Prot4Mom = trackP->Get4Momentum();
805 G4LorentzVector ProtSPposition = G4LorentzVector(trackP->GetPosition(), trackP->GetFormationTime());
806
807 for ( size_t j = 0; j < tracks->size(); ++j ) { // search for neutron
808
809 G4KineticTrack* trackN = (*tracks)[j];
810 if (! trackN ) continue;
811 if (trackN->GetDefinition() != neutron) continue;
812
813 G4LorentzVector Neut4Mom = trackN->Get4Momentum();
814 G4LorentzVector NeutSPposition = G4LorentzVector( trackN->GetPosition(), trackN->GetFormationTime()*hbarc/fermi);
815 G4double EffMass = (Prot4Mom + Neut4Mom).mag();
816
817 if ( EffMass <= MassCut ) { // && (EffDistance <= SpaceCut)) { // Create deuteron
818 G4KineticTrack* aDeuteron =
820 (trackP->GetFormationTime() + trackN->GetFormationTime())/2.0,
821 (trackP->GetPosition() + trackN->GetPosition() )/2.0,
822 ( Prot4Mom + Neut4Mom ));
823 aDeuteron->SetCreatorModelID(secID);
824 tracks->push_back(aDeuteron);
825 delete trackP; delete trackN;
826 (*tracks)[i] = nullptr; (*tracks)[j] = nullptr;
827 break;
828 }
829 }
830 }
831
832 // Find and remove null pointers created by decays above
833 for ( int jj = tracks->size()-1; jj >= 0; --jj ) {
834 if ( ! (*tracks)[jj] ) tracks->erase(tracks->begin()+jj);
835 }
836
837}
CLHEP::HepLorentzVector G4LorentzVector
static constexpr double fermi
Definition: G4SIunits.hh:83
double G4double
Definition: G4Types.hh:83
G4double GetFormationTime() const
void SetCreatorModelID(G4int id)
const G4ThreeVector & GetPosition() const
const G4ParticleDefinition * GetDefinition() const
const G4LorentzVector & Get4Momentum() const
float hbarc
Definition: hepunit.py:264

References DeltaM, deuteron, fermi, G4KineticTrack::Get4Momentum(), G4KineticTrack::GetDefinition(), G4KineticTrack::GetFormationTime(), G4ParticleDefinition::GetPDGMass(), G4KineticTrack::GetPosition(), source.hepunit::hbarc, neutron, proton, secID, and G4KineticTrack::SetCreatorModelID().

Referenced by PropagateNuclNucl().

◆ ModelDescription()

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

Reimplemented from G4HadronicInteraction.

Reimplemented in G4BinaryCascade, G4CascadeInterface, and G4INCLXXInterface.

Definition at line 50 of file G4VIntraNuclearTransportModel.cc.

51{
52 outFile << "G4VIntraNuclearTransportModel is abstract class.\n";
53 G4Exception("G4VIntraNuclearTransportModel::ModelDescription()","G4VINT01",
55 "G4VIntraNuclearTransportModel is abstract class, no description available");
56}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35

References FatalException, and G4Exception().

◆ 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 G4GeneratorPrecompoundInterface::operator!= ( G4GeneratorPrecompoundInterface right)
inlineprivate

Definition at line 90 of file G4GeneratorPrecompoundInterface.hh.

90{return (this != &right);}

◆ operator=()

const G4GeneratorPrecompoundInterface & G4GeneratorPrecompoundInterface::operator= ( const G4GeneratorPrecompoundInterface right)
private

◆ 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 G4GeneratorPrecompoundInterface::operator== ( G4GeneratorPrecompoundInterface right)
inlineprivate

Definition at line 89 of file G4GeneratorPrecompoundInterface.hh.

89{return (this == &right);}

◆ Propagate()

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

Implements G4VIntraNuclearTransportModel.

Definition at line 133 of file G4GeneratorPrecompoundInterface.cc.

135{
136 #ifdef debugPrecoInt
137 G4cout<<G4endl<<"G4GeneratorPrecompoundInterface::Propagate"<<G4endl;
138 G4cout<<"Target A and Z "<<theNucleus->GetMassNumber()<<" "<<theNucleus->GetCharge()<<G4endl;
139 G4cout<<"Directly produced particles number "<<theSecondaries->size()<<G4endl;
140 #endif
141
143
144 // decay the strong resonances
145 G4DecayKineticTracks decay(theSecondaries);
146 #ifdef debugPrecoInt
147 G4cout<<"Final stable particles number "<<theSecondaries->size()<<G4endl;
148 #endif
149
150 // prepare the fragment (it is assumed that target nuclei are never hypernuclei)
151 G4int anA=theNucleus->GetMassNumber();
152 G4int aZ=theNucleus->GetCharge();
153// G4double TargetNucleusMass = G4NucleiProperties::GetNuclearMass(anA, aZ);
154
155 G4int numberOfEx = 0;
156 G4int numberOfCh = 0;
157 G4int numberOfHoles = 0;
158
159 G4double R = theNucleus->GetNuclearRadius();
160
161 G4LorentzVector captured4Momentum(0.,0.,0.,0.);
162 G4LorentzVector Residual4Momentum(0.,0.,0.,0.); // TargetNucleusMass is not need at the moment
163 G4LorentzVector Secondary4Momentum(0.,0.,0.,0.);
164
165 // loop over secondaries
166 G4KineticTrackVector::iterator iter;
167 for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
168 {
169 const G4ParticleDefinition* part = (*iter)->GetDefinition();
170 G4double e = (*iter)->Get4Momentum().e();
171 G4double mass = (*iter)->Get4Momentum().mag();
172 G4ThreeVector mom = (*iter)->Get4Momentum().vect();
173 if((part != proton && part != neutron) ||
174 ((*iter)->GetPosition().mag() > R)) {
175 G4ReactionProduct * theNew = new G4ReactionProduct(part);
176 theNew->SetMomentum(mom);
177 theNew->SetTotalEnergy(e);
178 theNew->SetCreatorModelID((*iter)->GetCreatorModelID());
179 theTotalResult->push_back(theNew);
180 Secondary4Momentum += (*iter)->Get4Momentum();
181 #ifdef debugPrecoInt
182 G4cout<<"Secondary 4Mom "<<part->GetParticleName()<<" "<<(*iter)->Get4Momentum()<<" "
183 <<(*iter)->Get4Momentum().mag()<<G4endl;
184 #endif
185 } else {
186 if( e-mass > -CaptureThreshold*G4Log( G4UniformRand()) ) {
187 G4ReactionProduct * theNew = new G4ReactionProduct(part);
188 theNew->SetMomentum(mom);
189 theNew->SetTotalEnergy(e);
190 theNew->SetCreatorModelID((*iter)->GetCreatorModelID());
191 theTotalResult->push_back(theNew);
192 Secondary4Momentum += (*iter)->Get4Momentum();
193 #ifdef debugPrecoInt
194 G4cout<<"Secondary 4Mom "<<part->GetParticleName()<<" "<<(*iter)->Get4Momentum()<<" "
195 <<(*iter)->Get4Momentum().mag()<<G4endl;
196 #endif
197 } else {
198 // within the nucleus, neutron or proton
199 // now calculate A, Z of the fragment, momentum, number of exciton states
200 ++anA;
201 ++numberOfEx;
202 G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
203 aZ += Z;
204 numberOfCh += Z;
205 captured4Momentum += (*iter)->Get4Momentum();
206 #ifdef debugPrecoInt
207 G4cout<<"Captured 4Mom "<<part->GetParticleName()<<(*iter)->Get4Momentum()<<G4endl;
208 #endif
209 }
210 }
211 delete (*iter);
212 }
213 delete theSecondaries;
214
215 // loop over wounded nucleus
216 G4Nucleon * theCurrentNucleon =
217 theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
218 while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
219 {
220 if(theCurrentNucleon->AreYouHit()) {
221 ++numberOfHoles;
222 ++numberOfEx;
223 --anA;
224 aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/eplus + 0.1);
225
226 Residual4Momentum -= theCurrentNucleon->Get4Momentum();
227 }
228 theCurrentNucleon = theNucleus->GetNextNucleon();
229 }
230
231 #ifdef debugPrecoInt
232 G4cout<<G4endl;
233 G4cout<<"Secondary 4Mom "<<Secondary4Momentum<<G4endl;
234 G4cout<<"Captured 4Mom "<<captured4Momentum<<G4endl;
235 G4cout<<"Sec + Captured "<<Secondary4Momentum+captured4Momentum<<G4endl;
236 G4cout<<"Residual4Mom "<<Residual4Momentum<<G4endl;
237 G4cout<<"Sum 4 momenta "
238 <<Secondary4Momentum + captured4Momentum + Residual4Momentum <<G4endl;
239 #endif
240
241 // Check that we use QGS model; loop over wounded nucleus
242 G4bool QGSM(false);
243 theCurrentNucleon = theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
244 while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
245 {
246 if(theCurrentNucleon->AreYouHit())
247 {
248 if(theCurrentNucleon->Get4Momentum().mag() <
249 theCurrentNucleon->GetDefinition()->GetPDGMass()) QGSM=true;
250 }
251 theCurrentNucleon = theNucleus->GetNextNucleon();
252 }
253
254 #ifdef debugPrecoInt
255 if(!QGSM){
256 G4cout<<G4endl;
257 G4cout<<"Residual A and Z "<<anA<<" "<<aZ<<G4endl;
258 G4cout<<"Residual 4Mom "<<Residual4Momentum<<G4endl;
259 if(numberOfEx == 0)
260 {G4cout<<"Residual 4Mom = 0 means that there were not wounded and captured nucleons"<<G4endl;}
261 }
262 #endif
263
264 if(anA == 0) return theTotalResult;
265
266 G4LorentzVector exciton4Momentum(0.,0.,0.,0.);
267 if(anA >= aZ)
268 {
269 if(!QGSM)
270 { // FTF model was used
272
273// G4LorentzVector exciton4Momentum = Residual4Momentum + captured4Momentum;
274 exciton4Momentum = Residual4Momentum + captured4Momentum;
275//exciton4Momentum.setE(std::sqrt(exciton4Momentum.vect().mag2()+sqr(fMass)));
276 G4double ActualMass = exciton4Momentum.mag();
277 if(ActualMass <= fMass ) {
278 exciton4Momentum.setE(std::sqrt(exciton4Momentum.vect().mag2()+sqr(fMass)));
279 }
280
281 #ifdef debugPrecoInt
282 G4double exEnergy = 0.0;
283 if(ActualMass <= fMass ) {exEnergy = 0.;}
284 else {exEnergy = ActualMass - fMass;}
285 G4cout<<"Ground state residual Mass "<<fMass<<" E* "<<exEnergy<<G4endl;
286 #endif
287 }
288 else
289 { // QGS model was used
290 G4double InitialTargetMass =
291 G4NucleiProperties::GetNuclearMass(theNucleus->GetMassNumber(), theNucleus->GetCharge());
292
293 exciton4Momentum =
294 GetPrimaryProjectile()->Get4Momentum() + G4LorentzVector(0.,0.,0.,InitialTargetMass)
295 -Secondary4Momentum;
296
298 G4double ActualMass = exciton4Momentum.mag();
299
300 #ifdef debugPrecoInt
301 G4cout<<G4endl;
302 G4cout<<"Residual A and Z "<<anA<<" "<<aZ<<G4endl;
303 G4cout<<"Residual4Momentum "<<exciton4Momentum<<G4endl;
304 G4cout<<"ResidualMass, GroundStateMass and E* "<<ActualMass<<" "<<fMass<<" "
305 <<ActualMass - fMass<<G4endl;
306 #endif
307
308 if(ActualMass - fMass < 0.)
309 {
310 G4double ResE = std::sqrt(exciton4Momentum.vect().mag2() + sqr(fMass+10*MeV));
311 exciton4Momentum.setE(ResE);
312 #ifdef debugPrecoInt
313 G4cout<<"ActualMass - fMass < 0. "<<ActualMass<<" "<<fMass<<" "<<ActualMass - fMass<<G4endl;
314 G4int Uzhi; G4cin>>Uzhi;
315 #endif
316 }
317 }
318
319 // Need to de-excite the remnant nucleus only if excitation energy > 0.
320 G4Fragment anInitialState(anA, aZ, exciton4Momentum);
321 anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
322 anInitialState.SetNumberOfCharged(numberOfCh);
323 anInitialState.SetNumberOfHoles(numberOfHoles);
324 anInitialState.SetCreatorModelID(secID);
325
326 G4ReactionProductVector * aPrecoResult =
327 theDeExcitation->DeExcite(anInitialState);
328 // fill pre-compound part into the result, and return
329 #ifdef debugPrecoInt
330 G4cout<<"Target fragment number "<<aPrecoResult->size()<<G4endl;
331 #endif
332 for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
333 {
334 theTotalResult->push_back(aPrecoResult->operator[](ll));
335 #ifdef debugPrecoInt
336 G4cout<<"Fragment "<<ll<<" "
337 <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
338 <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
339 <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
340 <<aPrecoResult->operator[](ll)->GetDefinition()->GetPDGMass()<<G4endl;
341 #endif
342 }
343 delete aPrecoResult;
344 }
345
346 return theTotalResult;
347}
G4double G4Log(G4double x)
Definition: G4Log.hh:226
std::vector< G4ReactionProduct * > G4ReactionProductVector
static constexpr double eplus
Definition: G4SIunits.hh:184
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
#define G4cin
Definition: G4ios.hh:56
#define G4UniformRand()
Definition: Randomize.hh:52
const G4LorentzVector & Get4Momentum() const
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4bool AreYouHit() const
Definition: G4Nucleon.hh:98
virtual const G4LorentzVector & Get4Momentum() const
Definition: G4Nucleon.hh:72
virtual const G4ParticleDefinition * GetDefinition() const
Definition: G4Nucleon.hh:86
G4double GetPDGCharge() const
const G4String & GetParticleName() const
void SetMomentum(const G4double x, const G4double y, const G4double z)
void SetTotalEnergy(const G4double en)
void SetCreatorModelID(const G4int mod)
virtual G4Nucleon * GetNextNucleon()=0
virtual G4int GetCharge()=0
virtual G4bool StartLoop()=0
virtual G4double GetNuclearRadius()=0
virtual G4int GetMassNumber()=0
const G4HadProjectile * GetPrimaryProjectile() const
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)=0
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
T sqr(const T &x)
Definition: templates.hh:128

References G4Nucleon::AreYouHit(), CaptureThreshold, G4INCL::ClusterDecay::decay(), G4VPreCompoundModel::DeExcite(), eplus, G4cin, G4cout, G4endl, G4Log(), G4UniformRand, G4HadProjectile::Get4Momentum(), G4Nucleon::Get4Momentum(), G4V3DNucleus::GetCharge(), G4Nucleon::GetDefinition(), G4V3DNucleus::GetMassNumber(), G4V3DNucleus::GetNextNucleon(), G4NucleiProperties::GetNuclearMass(), G4V3DNucleus::GetNuclearRadius(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGCharge(), G4ParticleDefinition::GetPDGMass(), G4VIntraNuclearTransportModel::GetPrimaryProjectile(), CLHEP::HepLorentzVector::mag(), CLHEP::Hep3Vector::mag2(), MeV, neutron, proton, secID, G4ReactionProduct::SetCreatorModelID(), G4Fragment::SetCreatorModelID(), CLHEP::HepLorentzVector::setE(), G4ReactionProduct::SetMomentum(), G4Fragment::SetNumberOfCharged(), G4Fragment::SetNumberOfHoles(), G4Fragment::SetNumberOfParticles(), G4ReactionProduct::SetTotalEnergy(), sqr(), G4V3DNucleus::StartLoop(), G4VIntraNuclearTransportModel::theDeExcitation, CLHEP::HepLorentzVector::vect(), and Z.

◆ PropagateModelDescription()

void G4GeneratorPrecompoundInterface::PropagateModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 359 of file G4GeneratorPrecompoundInterface.cc.

360{
361 outFile << "G4GeneratorPrecompoundInterface interfaces a high\n"
362 << "energy model through the wounded nucleus to precompound de-excitation.\n"
363 << "Low energy protons and neutron present among secondaries produced by \n"
364 << "the high energy generator and within the nucleus are captured. The wounded\n"
365 << "nucleus and the captured particles form an excited nuclear fragment. This\n"
366 << "fragment is passed to the Geant4 pre-compound model for de-excitation.\n"
367 << "Nuclear de-excitation:\n";
368 // preco
369
370}

◆ PropagateNuclNucl()

G4ReactionProductVector * G4GeneratorPrecompoundInterface::PropagateNuclNucl ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus,
G4V3DNucleus theProjectileNucleus 
)
virtual

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 373 of file G4GeneratorPrecompoundInterface.cc.

376{
377#ifdef debugPrecoInt
378 G4cout<<G4endl<<"G4GeneratorPrecompoundInterface::PropagateNuclNucl "<<G4endl;
379 G4cout<<"Projectile A and Z (and numberOfLambdas) "<<theProjectileNucleus->GetMassNumber()<<" "
380 <<theProjectileNucleus->GetCharge()<<" ("
381 <<theProjectileNucleus->GetNumberOfLambdas()<<")"<<G4endl;
382 G4cout<<"Target A and Z "<<theNucleus->GetMassNumber()<<" "
383 <<theNucleus->GetCharge()<<" ("
384 <<theNucleus->GetNumberOfLambdas()<<")"<<G4endl;
385 G4cout<<"Directly produced particles number "<<theSecondaries->size()<<G4endl;
386 G4cout<<"Projectile 4Mom and mass "<<GetPrimaryProjectile()->Get4Momentum()<<" "
388#endif
389
390 // prepare the target residual (assumed to be never a hypernucleus)
391 G4int anA=theNucleus->GetMassNumber();
392 G4int aZ=theNucleus->GetCharge();
393 //G4int aL=theNucleus->GetNumberOfLambdas(); // Should be 0
394 G4int numberOfEx = 0;
395 G4int numberOfCh = 0;
396 G4int numberOfHoles = 0;
397 G4double exEnergy = 0.0;
398 G4double R = theNucleus->GetNuclearRadius();
399 G4LorentzVector Target4Momentum(0.,0.,0.,0.);
400
401 // loop over the wounded target nucleus
402 G4Nucleon * theCurrentNucleon =
403 theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
404 while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
405 {
406 if(theCurrentNucleon->AreYouHit()) {
407 ++numberOfHoles;
408 ++numberOfEx;
409 --anA;
410 aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
411 eplus + 0.1);
412 exEnergy += theCurrentNucleon->GetBindingEnergy();
413 Target4Momentum -=theCurrentNucleon->Get4Momentum();
414 }
415 theCurrentNucleon = theNucleus->GetNextNucleon();
416 }
417
418#ifdef debugPrecoInt
419 G4cout<<"Residual Target A Z (numberOfLambdas) E* 4mom "<<anA<<" "<<aZ<<" (0"//<<aL
420 <<") "<<exEnergy<<" "<<Target4Momentum<<G4endl;
421#endif
422
423 // prepare the projectile residual - which can be a hypernucleus or anti-hypernucleus
424
425 G4bool ProjectileIsAntiNucleus=
427
429
430 G4int anAb=theProjectileNucleus->GetMassNumber();
431 G4int aZb=theProjectileNucleus->GetCharge();
432 G4int aLb=theProjectileNucleus->GetNumberOfLambdas(); // Non negative number of (anti-)lambdas in (anti-)nucleus
433 G4int numberOfExB = 0;
434 G4int numberOfChB = 0;
435 G4int numberOfHolesB = 0;
436 G4double exEnergyB = 0.0;
437 G4double Rb = theProjectileNucleus->GetNuclearRadius();
438 G4LorentzVector Projectile4Momentum(0.,0.,0.,0.);
439
440 // loop over the wounded projectile nucleus or anti-nucleus
441 theCurrentNucleon =
442 theProjectileNucleus->StartLoop() ? theProjectileNucleus->GetNextNucleon() : 0;
443 while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
444 {
445 if(theCurrentNucleon->AreYouHit()) {
446 ++numberOfHolesB;
447 ++numberOfExB;
448 --anAb;
449 if(!ProjectileIsAntiNucleus) {
450 aZb -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
451 eplus + 0.1);
452 if (theCurrentNucleon->GetParticleType()==G4Lambda::Definition()) --aLb;
453 } else {
454 aZb += G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
455 eplus - 0.1);
456 if (theCurrentNucleon->GetParticleType()==G4AntiLambda::Definition()) --aLb;
457 }
458 exEnergyB += theCurrentNucleon->GetBindingEnergy();
459 Projectile4Momentum -=theCurrentNucleon->Get4Momentum();
460 }
461 theCurrentNucleon = theProjectileNucleus->GetNextNucleon();
462 }
463
464 G4bool ExistTargetRemnant = G4double (numberOfHoles) <
465 0.3* G4double (numberOfHoles + anA);
466 G4bool ExistProjectileRemnant= G4double (numberOfHolesB) <
467 0.3*G4double (numberOfHolesB + anAb);
468
469#ifdef debugPrecoInt
470 G4cout<<"Projectile residual A Z (numberOfLambdas) E* 4mom "<<anAb<<" "<<aZb<<" ("<<aLb
471 <<") "<<exEnergyB<<" "<<Projectile4Momentum<<G4endl;
472 G4cout<<" ExistTargetRemnant ExistProjectileRemnant "
473 <<ExistTargetRemnant<<" "<< ExistProjectileRemnant<<G4endl;
474#endif
475 //-----------------------------------------------------------------------------
476 // decay the strong resonances
478 G4DecayKineticTracks decay(theSecondaries);
479
480 MakeCoalescence(theSecondaries);
481
482 #ifdef debugPrecoInt
483 G4cout<<"Secondary stable particles number "<<theSecondaries->size()<<G4endl;
484 #endif
485
486#ifdef debugPrecoInt
487 G4LorentzVector secondary4Momemtum(0,0,0,0);
488 G4int SecondrNum(0);
489#endif
490
491 // Loop over secondaries.
492 // We are assuming that only protons and neutrons - for nuclei -
493 // and only antiprotons and antineutrons - for antinuclei - can be absorbed,
494 // not instead lambdas (or hyperons more generally) - for nuclei - or anti-lambdas
495 // (or anti-hyperons more generally) - for antinuclei. This is a simplification,
496 // to be eventually reviewed later on, in particular when generic hypernuclei and
497 // anti-hypernuclei are introduced, instead of the few light hypernuclei and
498 // anti-hypernuclei which currently exist in Geant4.
499 G4KineticTrackVector::iterator iter;
500 for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
501 {
502 const G4ParticleDefinition* part = (*iter)->GetDefinition();
503 G4LorentzVector aTrack4Momentum=(*iter)->Get4Momentum();
504
505 if( part != proton && part != neutron &&
506 (part != ANTIproton && ProjectileIsAntiNucleus) &&
507 (part != ANTIneutron && ProjectileIsAntiNucleus) )
508 {
509 G4ReactionProduct * theNew = new G4ReactionProduct(part);
510 theNew->SetMomentum(aTrack4Momentum.vect());
511 theNew->SetTotalEnergy(aTrack4Momentum.e());
512 theNew->SetCreatorModelID((*iter)->GetCreatorModelID());
513 theTotalResult->push_back(theNew);
514#ifdef debugPrecoInt
515 SecondrNum++;
516 secondary4Momemtum += (*iter)->Get4Momentum();
517 G4cout<<"Secondary "<<SecondrNum<<" "
518 <<theNew->GetDefinition()->GetParticleName()<<" "
519 <<theNew->GetMomentum()<<" "<<theNew->GetTotalEnergy()<<G4endl;
520
521#endif
522 delete (*iter);
523 continue;
524 }
525
526 G4bool CanBeCapturedByTarget = false;
527 if( part == proton || part == neutron)
528 {
529 CanBeCapturedByTarget = ExistTargetRemnant &&
531 (aTrack4Momentum + Target4Momentum).mag() -
532 aTrack4Momentum.mag() - Target4Momentum.mag()) &&
533 ((*iter)->GetPosition().mag() < R);
534 }
535 // ---------------------------
536 G4LorentzVector Position((*iter)->GetPosition(), (*iter)->GetFormationTime());
537 Position.boost(bst);
538
539 G4bool CanBeCapturedByProjectile = false;
540
541 if( !ProjectileIsAntiNucleus &&
542 ( part == proton || part == neutron))
543 {
544 CanBeCapturedByProjectile = ExistProjectileRemnant &&
546 (aTrack4Momentum + Projectile4Momentum).mag() -
547 aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
548 (Position.vect().mag() < Rb);
549 }
550
551 if( ProjectileIsAntiNucleus &&
552 ( part == ANTIproton || part == ANTIneutron))
553 {
554 CanBeCapturedByProjectile = ExistProjectileRemnant &&
556 (aTrack4Momentum + Projectile4Momentum).mag() -
557 aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
558 (Position.vect().mag() < Rb);
559 }
560
561 if(CanBeCapturedByTarget && CanBeCapturedByProjectile)
562 {
563 if(G4UniformRand() < 0.5)
564 { CanBeCapturedByTarget = true; CanBeCapturedByProjectile = false;}
565 else
566 { CanBeCapturedByTarget = false; CanBeCapturedByProjectile = true;}
567 }
568
569 if(CanBeCapturedByTarget)
570 {
571 // within the target nucleus, neutron or proton
572 // now calculate A, Z of the fragment, momentum,
573 // number of exciton states
574#ifdef debugPrecoInt
575 G4cout<<"Track is CapturedByTarget "<<" "<<part->GetParticleName()<<" "
576 <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
577#endif
578 ++anA;
579 ++numberOfEx;
580 G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
581 aZ += Z;
582 numberOfCh += Z;
583 Target4Momentum +=aTrack4Momentum;
584 delete (*iter);
585 } else if(CanBeCapturedByProjectile)
586 {
587 // within the projectile nucleus, neutron or proton
588 // now calculate A, Z of the fragment, momentum,
589 // number of exciton states
590#ifdef debugPrecoInt
591 G4cout<<"Track is CapturedByProjectile"<<" "<<part->GetParticleName()<<" "
592 <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
593#endif
594 ++anAb;
595 ++numberOfExB;
596 G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
597 if( ProjectileIsAntiNucleus ) Z=-Z;
598 aZb += Z;
599 numberOfChB += Z;
600 Projectile4Momentum +=aTrack4Momentum;
601 delete (*iter);
602 } else
603 { // the track is not captured
604 G4ReactionProduct * theNew = new G4ReactionProduct(part);
605 theNew->SetMomentum(aTrack4Momentum.vect());
606 theNew->SetTotalEnergy(aTrack4Momentum.e());
607 theNew->SetCreatorModelID((*iter)->GetCreatorModelID());
608 theTotalResult->push_back(theNew);
609
610#ifdef debugPrecoInt
611 SecondrNum++;
612 secondary4Momemtum += (*iter)->Get4Momentum();
613/*
614 G4cout<<"Secondary "<<SecondrNum<<" "
615 <<theNew->GetDefinition()->GetParticleName()<<" "
616 <<secondary4Momemtum<<G4endl;
617*/
618#endif
619 delete (*iter);
620 continue;
621 }
622 }
623 delete theSecondaries;
624 //-----------------------------------------------------
625
626 #ifdef debugPrecoInt
627 G4cout<<"Final target residual A Z (numberOfLambdas) E* 4mom "<<anA<<" "<<aZ<<" (0"//<<aL
628 <<") "<<exEnergy<<" "<<Target4Momentum<<G4endl;
629 #endif
630
631 if(0!=anA )
632 {
633 // We assume that the target residual is never a hypernucleus
635
636 if((anA == theNucleus->GetMassNumber()) && (exEnergy <= 0.))
637 {Target4Momentum.setE(fMass);}
638
639 G4double RemnMass=Target4Momentum.mag();
640
641 if(RemnMass < fMass)
642 {
643 RemnMass=fMass + exEnergy;
644 Target4Momentum.setE(std::sqrt(Target4Momentum.vect().mag2() +
645 RemnMass*RemnMass));
646 } else
647 { exEnergy=RemnMass-fMass;}
648
649 if( exEnergy < 0.) exEnergy=0.;
650
651 // Need to de-excite the remnant nucleus
652 G4Fragment anInitialState(anA, aZ, Target4Momentum);
653 anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
654 anInitialState.SetNumberOfCharged(numberOfCh);
655 anInitialState.SetNumberOfHoles(numberOfHoles);
656 anInitialState.SetCreatorModelID(secID);
657
658 G4ReactionProductVector * aPrecoResult =
659 theDeExcitation->DeExcite(anInitialState);
660
661 #ifdef debugPrecoInt
662 G4cout<<"Target fragment number "<<aPrecoResult->size()<<G4endl;
663 #endif
664
665 // fill pre-compound part into the result, and return
666 for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
667 {
668 theTotalResult->push_back(aPrecoResult->operator[](ll));
669 #ifdef debugPrecoInt
670 G4cout<<"Target fragment "<<ll<<" "
671 <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
672 <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
673 <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
674 <<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
675 #endif
676 }
677 delete aPrecoResult;
678 }
679
680 //-----------------------------------------------------
681 if((anAb == theProjectileNucleus->GetMassNumber())&& (exEnergyB <= 0.))
682 {Projectile4Momentum = GetPrimaryProjectile()->Get4Momentum();}
683
684 #ifdef debugPrecoInt
685 G4cout<<"Final projectile residual A Z (numberOfLambdas) E* Pmom Pmag2 "<<anAb<<" "<<aZb<<" ("
686 <<aLb<<") "<<exEnergyB<<" "<<Projectile4Momentum<<" "
687 <<Projectile4Momentum.mag2()<<G4endl;
688 #endif
689
690 if(0!=anAb)
691 {
692 // The projectile residual can be a hypernucleus or anti-hypernucleus
693 G4double fMass = 0.0;
694 if ( aLb > 0 ) {
695 fMass = G4HyperNucleiProperties::GetNuclearMass(anAb, aZb, aLb);
696 } else {
697 fMass = G4NucleiProperties::GetNuclearMass(anAb, aZb);
698 }
699 G4double RemnMass=Projectile4Momentum.mag();
700
701 if(RemnMass < fMass)
702 {
703 RemnMass=fMass + exEnergyB;
704 Projectile4Momentum.setE(std::sqrt(Projectile4Momentum.vect().mag2() +
705 RemnMass*RemnMass));
706 } else
707 { exEnergyB=RemnMass-fMass;}
708
709 if( exEnergyB < 0.) exEnergyB=0.;
710
711 G4ThreeVector bstToCM =Projectile4Momentum.findBoostToCM();
712 Projectile4Momentum.boost(bstToCM);
713
714 // Need to de-excite the remnant nucleus
715 G4Fragment anInitialState(anAb, aZb, aLb, Projectile4Momentum);
716 anInitialState.SetNumberOfParticles(numberOfExB-numberOfHolesB);
717 anInitialState.SetNumberOfCharged(numberOfChB);
718 anInitialState.SetNumberOfHoles(numberOfHolesB);
719 anInitialState.SetCreatorModelID(secID);
720
721 G4ReactionProductVector * aPrecoResult =
722 theDeExcitation->DeExcite(anInitialState);
723
724 #ifdef debugPrecoInt
725 G4cout<<"Projectile fragment number "<<aPrecoResult->size()<<G4endl;
726 #endif
727
728 // fill pre-compound part into the result, and return
729 for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
730 {
731 G4LorentzVector tmp=G4LorentzVector(aPrecoResult->operator[](ll)->GetMomentum(),
732 aPrecoResult->operator[](ll)->GetTotalEnergy());
733 tmp.boost(-bstToCM); // Transformation to the system of original remnant
734 aPrecoResult->operator[](ll)->SetMomentum(tmp.vect());
735 aPrecoResult->operator[](ll)->SetTotalEnergy(tmp.e());
736
737 if(ProjectileIsAntiNucleus)
738 {
739 const G4ParticleDefinition * aFragment=aPrecoResult->operator[](ll)->GetDefinition();
740 const G4ParticleDefinition * LastFragment=aFragment;
741 if (aFragment == proton) {LastFragment=G4AntiProton::AntiProtonDefinition();}
742 else if(aFragment == neutron) {LastFragment=G4AntiNeutron::AntiNeutronDefinition();}
743 else if(aFragment == deuteron){LastFragment=G4AntiDeuteron::AntiDeuteronDefinition();}
744 else if(aFragment == triton) {LastFragment=G4AntiTriton::AntiTritonDefinition();}
745 else if(aFragment == He3) {LastFragment=G4AntiHe3::AntiHe3Definition();}
746 else if(aFragment == He4) {LastFragment=G4AntiAlpha::AntiAlphaDefinition();}
747 else {}
748
749 if (aLb > 0) { // Anti-hypernucleus
750 if (aFragment == G4HyperTriton::Definition()) {
751 LastFragment=G4AntiHyperTriton::Definition();
752 } else if (aFragment == G4HyperH4::Definition()) {
753 LastFragment=G4AntiHyperH4::Definition();
754 } else if (aFragment == G4HyperAlpha::Definition()) {
755 LastFragment=G4AntiHyperAlpha::Definition();
756 } else if (aFragment == G4HyperHe5::Definition()) {
757 LastFragment=G4AntiHyperHe5::Definition();
758 } else if (aFragment == G4DoubleHyperH4::Definition()) {
759 LastFragment=G4AntiDoubleHyperH4::Definition();
760 } else if (aFragment == G4DoubleHyperDoubleNeutron::Definition()) {
762 }
763 }
764
765 aPrecoResult->operator[](ll)->SetDefinitionAndUpdateE(LastFragment);
766 }
767
768 #ifdef debugPrecoInt
769 G4cout<<"Projectile fragment "<<ll<<" "
770 <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
771 <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
772 <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
773 <<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
774 #endif
775
776 theTotalResult->push_back(aPrecoResult->operator[](ll));
777 }
778
779 delete aPrecoResult;
780 }
781
782 return theTotalResult;
783}
Hep3Vector boostVector() const
HepLorentzVector & boost(double, double, double)
Hep3Vector vect() const
static G4AntiAlpha * AntiAlphaDefinition()
Definition: G4AntiAlpha.cc:83
static G4AntiDeuteron * AntiDeuteronDefinition()
static G4AntiDoubleHyperDoubleNeutron * Definition()
static G4AntiDoubleHyperH4 * Definition()
static G4AntiHe3 * AntiHe3Definition()
Definition: G4AntiHe3.cc:88
static G4AntiHyperAlpha * Definition()
static G4AntiHyperH4 * Definition()
static G4AntiHyperHe5 * Definition()
static G4AntiHyperTriton * Definition()
static G4AntiLambda * Definition()
Definition: G4AntiLambda.cc:52
static G4AntiNeutron * AntiNeutronDefinition()
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:87
static G4AntiTriton * AntiTritonDefinition()
Definition: G4AntiTriton.cc:88
static G4DoubleHyperDoubleNeutron * Definition()
static G4DoubleHyperH4 * Definition()
void MakeCoalescence(G4KineticTrackVector *theSecondaries)
const G4ParticleDefinition * GetDefinition() const
static G4HyperAlpha * Definition()
Definition: G4HyperAlpha.cc:45
static G4HyperH4 * Definition()
Definition: G4HyperH4.cc:45
static G4HyperHe5 * Definition()
Definition: G4HyperHe5.cc:45
static G4double GetNuclearMass(G4int A, G4int Z, G4int L)
static G4HyperTriton * Definition()
static G4Lambda * Definition()
Definition: G4Lambda.cc:52
const G4ParticleDefinition * GetParticleType() const
Definition: G4Nucleon.hh:85
G4double GetBindingEnergy() const
Definition: G4Nucleon.hh:75
const G4ParticleDefinition * GetDefinition() const
G4double GetTotalEnergy() const
G4ThreeVector GetMomentum() const
virtual G4int GetNumberOfLambdas()=0

References G4AntiAlpha::AntiAlphaDefinition(), G4AntiDeuteron::AntiDeuteronDefinition(), G4AntiHe3::AntiHe3Definition(), ANTIneutron, G4AntiNeutron::AntiNeutronDefinition(), ANTIproton, G4AntiProton::AntiProtonDefinition(), G4AntiTriton::AntiTritonDefinition(), G4Nucleon::AreYouHit(), CLHEP::HepLorentzVector::boost(), CLHEP::HepLorentzVector::boostVector(), CaptureThreshold, G4INCL::ClusterDecay::decay(), G4VPreCompoundModel::DeExcite(), G4AntiLambda::Definition(), G4Lambda::Definition(), G4AntiDoubleHyperDoubleNeutron::Definition(), G4AntiDoubleHyperH4::Definition(), G4AntiHyperAlpha::Definition(), G4AntiHyperH4::Definition(), G4AntiHyperHe5::Definition(), G4AntiHyperTriton::Definition(), G4DoubleHyperDoubleNeutron::Definition(), G4DoubleHyperH4::Definition(), G4HyperAlpha::Definition(), G4HyperH4::Definition(), G4HyperHe5::Definition(), G4HyperTriton::Definition(), deuteron, CLHEP::HepLorentzVector::e(), eplus, CLHEP::HepLorentzVector::findBoostToCM(), G4cout, G4endl, G4Log(), G4UniformRand, G4HadProjectile::Get4Momentum(), G4Nucleon::Get4Momentum(), G4ParticleDefinition::GetBaryonNumber(), G4Nucleon::GetBindingEnergy(), G4V3DNucleus::GetCharge(), G4HadProjectile::GetDefinition(), G4Nucleon::GetDefinition(), G4ReactionProduct::GetDefinition(), G4V3DNucleus::GetMassNumber(), G4ReactionProduct::GetMomentum(), G4V3DNucleus::GetNextNucleon(), G4NucleiProperties::GetNuclearMass(), G4HyperNucleiProperties::GetNuclearMass(), G4V3DNucleus::GetNuclearRadius(), G4V3DNucleus::GetNumberOfLambdas(), G4ParticleDefinition::GetParticleName(), G4Nucleon::GetParticleType(), G4ParticleDefinition::GetPDGCharge(), G4VIntraNuclearTransportModel::GetPrimaryProjectile(), G4ReactionProduct::GetTotalEnergy(), He3, He4, CLHEP::HepLorentzVector::mag(), CLHEP::HepLorentzVector::mag2(), CLHEP::Hep3Vector::mag2(), MakeCoalescence(), neutron, proton, secID, G4ReactionProduct::SetCreatorModelID(), G4Fragment::SetCreatorModelID(), CLHEP::HepLorentzVector::setE(), G4ReactionProduct::SetMomentum(), G4Fragment::SetNumberOfCharged(), G4Fragment::SetNumberOfHoles(), G4Fragment::SetNumberOfParticles(), G4ReactionProduct::SetTotalEnergy(), G4V3DNucleus::StartLoop(), G4VIntraNuclearTransportModel::theDeExcitation, triton, CLHEP::HepLorentzVector::vect(), and Z.

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

◆ SetCaptureThreshold()

void G4GeneratorPrecompoundInterface::SetCaptureThreshold ( G4double  value)
inline

Definition at line 116 of file G4GeneratorPrecompoundInterface.hh.

117{
118 CaptureThreshold = value;
119}

References CaptureThreshold.

◆ SetDeExcitation()

void G4VIntraNuclearTransportModel::SetDeExcitation ( G4VPreCompoundModel ptr)
inlineinherited

◆ SetDeltaM()

void G4GeneratorPrecompoundInterface::SetDeltaM ( G4double  value)
inline

Definition at line 122 of file G4GeneratorPrecompoundInterface.hh.

123{
124 DeltaM = value;
125}

References DeltaM.

◆ SetDeltaR()

void G4GeneratorPrecompoundInterface::SetDeltaR ( G4double  value)
inline

Definition at line 128 of file G4GeneratorPrecompoundInterface.hh.

129{
130 DeltaR = value;
131}

References DeltaR.

◆ 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

Field Documentation

◆ ANTIdeuteron

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTIdeuteron
private

Definition at line 107 of file G4GeneratorPrecompoundInterface.hh.

Referenced by G4GeneratorPrecompoundInterface().

◆ ANTIHe3

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTIHe3
private

Definition at line 109 of file G4GeneratorPrecompoundInterface.hh.

Referenced by G4GeneratorPrecompoundInterface().

◆ ANTIHe4

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTIHe4
private

Definition at line 110 of file G4GeneratorPrecompoundInterface.hh.

Referenced by G4GeneratorPrecompoundInterface().

◆ ANTIneutron

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTIneutron
private

◆ ANTIproton

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTIproton
private

◆ ANTItriton

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::ANTItriton
private

Definition at line 108 of file G4GeneratorPrecompoundInterface.hh.

Referenced by G4GeneratorPrecompoundInterface().

◆ CaptureThreshold

G4double G4GeneratorPrecompoundInterface::CaptureThreshold
private

◆ DeltaM

G4double G4GeneratorPrecompoundInterface::DeltaM
private

Definition at line 93 of file G4GeneratorPrecompoundInterface.hh.

Referenced by MakeCoalescence(), and SetDeltaM().

◆ DeltaR

G4double G4GeneratorPrecompoundInterface::DeltaR
private

Definition at line 94 of file G4GeneratorPrecompoundInterface.hh.

Referenced by SetDeltaR().

◆ deuteron

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::deuteron
private

◆ epCheckLevels

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

◆ He3

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::He3
private

◆ He4

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::He4
private

◆ isBlocked

G4bool G4HadronicInteraction::isBlocked
protectedinherited

◆ neutron

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::neutron
private

◆ proton

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::proton
private

◆ recoilEnergyThreshold

G4double G4HadronicInteraction::recoilEnergyThreshold
privateinherited

◆ registry

G4HadronicInteractionRegistry* G4HadronicInteraction::registry
privateinherited

◆ secID

G4int G4GeneratorPrecompoundInterface::secID
private

◆ the3DNucleus

G4V3DNucleus* G4VIntraNuclearTransportModel::the3DNucleus
protectedinherited

◆ theBlockedList

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

◆ theBlockedListElements

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

◆ theDeExcitation

G4VPreCompoundModel* G4VIntraNuclearTransportModel::theDeExcitation
protectedinherited

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

◆ thePrimaryProjectile

const G4HadProjectile* G4VIntraNuclearTransportModel::thePrimaryProjectile
protectedinherited

◆ theTransportModelName

G4String G4VIntraNuclearTransportModel::theTransportModelName
protectedinherited

◆ triton

const G4ParticleDefinition* G4GeneratorPrecompoundInterface::triton
private

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