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

#include <G4HESigmaMinusInelastic.hh>

+ Inheritance diagram for G4HESigmaMinusInelastic:

Public Member Functions

 G4HESigmaMinusInelastic (const G4String &name="G4HESigmaMinusInelastic")
 
 ~G4HESigmaMinusInelastic ()
 
virtual void ModelDescription (std::ostream &) const
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4int GetNumberOfSecondaries ()
 
void FirstIntInCasSigmaMinus (G4bool &inElastic, const G4double availableEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, const G4double atomicWeight)
 
- Public Member Functions inherited from G4HEInelastic
 G4HEInelastic (const G4String &modelName="HEInelastic")
 
 ~G4HEInelastic ()
 
void SetMaxNumberOfSecondaries (const G4int maxnumber)
 
void SetVerboseLevel (const G4int level)
 
void ForceEnergyConservation (G4bool energyConservation)
 
G4bool EnergyConservation (void)
 
virtual const std::pair< G4double, G4doubleGetFatalEnergyCheckLevels () const
 
G4double Amin (G4double a, G4double b)
 
G4double Amax (G4double a, G4double b)
 
G4int Imin (G4int a, G4int b)
 
G4int Imax (G4int a, G4int b)
 
void FillParticleChange (G4HEVector pv[], G4int aVecLength)
 
G4double pmltpc (G4int np, G4int nm, G4int nz, G4int n, G4double b, G4double c)
 
G4int Factorial (G4int n)
 
G4double NuclearInelasticity (G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber)
 
G4double NuclearExcitation (G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber, G4double &excitationEnergyCascade, G4double &excitationEnergyEvaporation)
 
void HighEnergyCascading (G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void HighEnergyClusterProduction (G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void TuningOfHighEnergyCascading (G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void MediumEnergyCascading (G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void MediumEnergyClusterProduction (G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void QuasiElasticScattering (G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
 
void ElasticScattering (G4bool &successful, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, G4double atomicWeight, G4double atomicNumber)
 
G4int rtmi (G4double *x, G4double xli, G4double xri, G4double eps, G4int iend, G4double aa, G4double bb, G4double cc, G4double dd, G4double rr)
 
G4double fctcos (G4double t, G4double aa, G4double bb, G4double cc, G4double dd, G4double rr)
 
void StrangeParticlePairProduction (const G4double availableEnergy, const G4double centerOfMassEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle)
 
G4double NBodyPhaseSpace (const G4double totalEnergy, const G4bool constantCrossSection, G4HEVector pv[], G4int &vecLen)
 
G4double NBodyPhaseSpace (G4int npart, G4HEVector pv[], G4double wmax, G4double wfcn, G4int maxtrial, G4int ntrial)
 
G4double gpdk (G4double a, G4double b, G4double c)
 
void QuickSort (G4double arr[], const G4int lidx, const G4int ridx)
 
G4double Alam (G4double a, G4double b, G4double c)
 
G4double CalculatePhaseSpaceWeight (G4int npart)
 
G4double normal (void)
 
G4double GammaRand (G4double avalue)
 
G4double Erlang (G4int mvalue)
 
G4int Poisson (G4double x)
 
void SetParticles (void)
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &, G4Nucleus &)
 
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)
 
const G4HadronicInteractionGetMyPointer () const
 
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
 
G4bool operator== (const G4HadronicInteraction &right) const
 
G4bool operator!= (const G4HadronicInteraction &right) 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
 

Public Attributes

G4int vecLength
 
- Public Attributes inherited from G4HEInelastic
G4int verboseLevel
 
G4int MAXPART
 
G4bool conserveEnergy
 
G4HEVector PionPlus
 
G4HEVector PionZero
 
G4HEVector PionMinus
 
G4HEVector KaonPlus
 
G4HEVector KaonZero
 
G4HEVector AntiKaonZero
 
G4HEVector KaonMinus
 
G4HEVector KaonZeroShort
 
G4HEVector KaonZeroLong
 
G4HEVector Proton
 
G4HEVector AntiProton
 
G4HEVector Neutron
 
G4HEVector AntiNeutron
 
G4HEVector Lambda
 
G4HEVector AntiLambda
 
G4HEVector SigmaPlus
 
G4HEVector SigmaZero
 
G4HEVector SigmaMinus
 
G4HEVector AntiSigmaPlus
 
G4HEVector AntiSigmaZero
 
G4HEVector AntiSigmaMinus
 
G4HEVector XiZero
 
G4HEVector XiMinus
 
G4HEVector AntiXiZero
 
G4HEVector AntiXiMinus
 
G4HEVector OmegaMinus
 
G4HEVector AntiOmegaMinus
 
G4HEVector Deuteron
 
G4HEVector Triton
 
G4HEVector Alpha
 
G4HEVector Gamma
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 51 of file G4HESigmaMinusInelastic.hh.

Constructor & Destructor Documentation

◆ G4HESigmaMinusInelastic()

G4HESigmaMinusInelastic::G4HESigmaMinusInelastic ( const G4String name = "G4HESigmaMinusInelastic")

Definition at line 42 of file G4HESigmaMinusInelastic.cc.

43 : G4HEInelastic(name)
44{
45 vecLength = 0;
46 theMinEnergy = 20*GeV;
47 theMaxEnergy = 10*TeV;
48 MAXPART = 2048;
49 verboseLevel = 0;
50 G4cout << "WARNING: model G4HESigmaMinusInelastic is being deprecated and will\n"
51 << "disappear in Geant4 version 10.0" << G4endl;
52}
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout

◆ ~G4HESigmaMinusInelastic()

G4HESigmaMinusInelastic::~G4HESigmaMinusInelastic ( )
inline

Definition at line 56 of file G4HESigmaMinusInelastic.hh.

56{};

Member Function Documentation

◆ ApplyYourself()

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

Implements G4HadronicInteraction.

Definition at line 69 of file G4HESigmaMinusInelastic.cc.

71{
72 G4HEVector* pv = new G4HEVector[MAXPART];
73 const G4HadProjectile* aParticle = &aTrack;
74 const G4double A = targetNucleus.GetA_asInt();
75 const G4double Z = targetNucleus.GetZ_asInt();
76 G4HEVector incidentParticle(aParticle);
77
78 G4double atomicNumber = Z;
79 G4double atomicWeight = A;
80
81 G4int incidentCode = incidentParticle.getCode();
82 G4double incidentMass = incidentParticle.getMass();
83 G4double incidentTotalEnergy = incidentParticle.getEnergy();
84
85 // G4double incidentTotalMomentum = incidentParticle.getTotalMomentum();
86 // DHW 19 May 2011: variable set but not used
87
88 G4double incidentKineticEnergy = incidentTotalEnergy - incidentMass;
89
90 if (incidentKineticEnergy < 1.)
91 G4cout << "GHESigmaMinusInelastic: incident energy < 1 GeV" << G4endl;
92
93 if (verboseLevel > 1) {
94 G4cout << "G4HESigmaMinusInelastic::ApplyYourself" << G4endl;
95 G4cout << "incident particle " << incidentParticle.getName()
96 << "mass " << incidentMass
97 << "kinetic energy " << incidentKineticEnergy
98 << G4endl;
99 G4cout << "target material with (A,Z) = ("
100 << atomicWeight << "," << atomicNumber << ")" << G4endl;
101 }
102
103 G4double inelasticity = NuclearInelasticity(incidentKineticEnergy,
104 atomicWeight, atomicNumber);
105 if (verboseLevel > 1)
106 G4cout << "nuclear inelasticity = " << inelasticity << G4endl;
107
108 incidentKineticEnergy -= inelasticity;
109
110 G4double excitationEnergyGNP = 0.;
111 G4double excitationEnergyDTA = 0.;
112
113 G4double excitation = NuclearExcitation(incidentKineticEnergy,
114 atomicWeight, atomicNumber,
115 excitationEnergyGNP,
116 excitationEnergyDTA);
117 if (verboseLevel > 1)
118 G4cout << "nuclear excitation = " << excitation << excitationEnergyGNP
119 << excitationEnergyDTA << G4endl;
120
121 incidentKineticEnergy -= excitation;
122 incidentTotalEnergy = incidentKineticEnergy + incidentMass;
123 // incidentTotalMomentum = std::sqrt( (incidentTotalEnergy-incidentMass)
124 // *(incidentTotalEnergy+incidentMass));
125 // DHW 19 May 2011: variable set but not used
126
127 G4HEVector targetParticle;
128 if (G4UniformRand() < atomicNumber/atomicWeight) {
129 targetParticle.setDefinition("Proton");
130 } else {
131 targetParticle.setDefinition("Neutron");
132 }
133
134 G4double targetMass = targetParticle.getMass();
135 G4double centerOfMassEnergy = std::sqrt(incidentMass*incidentMass
136 + targetMass*targetMass
137 + 2.0*targetMass*incidentTotalEnergy);
138 G4double availableEnergy = centerOfMassEnergy - targetMass - incidentMass;
139
140 G4bool inElastic = true;
141 vecLength = 0;
142
143 if (verboseLevel > 1)
144 G4cout << "ApplyYourself: CallFirstIntInCascade for particle "
145 << incidentCode << G4endl;
146
147 G4bool successful = false;
148
149 FirstIntInCasSigmaMinus(inElastic, availableEnergy, pv, vecLength,
150 incidentParticle, targetParticle, atomicWeight);
151
152 if (verboseLevel > 1)
153 G4cout << "ApplyYourself::StrangeParticlePairProduction" << G4endl;
154
155 if ((vecLength > 0) && (availableEnergy > 1.))
156 StrangeParticlePairProduction(availableEnergy, centerOfMassEnergy,
157 pv, vecLength,
158 incidentParticle, targetParticle);
159
160 HighEnergyCascading(successful, pv, vecLength,
161 excitationEnergyGNP, excitationEnergyDTA,
162 incidentParticle, targetParticle,
163 atomicWeight, atomicNumber);
164 if (!successful)
166 excitationEnergyGNP, excitationEnergyDTA,
167 incidentParticle, targetParticle,
168 atomicWeight, atomicNumber);
169 if (!successful)
170 MediumEnergyCascading(successful, pv, vecLength,
171 excitationEnergyGNP, excitationEnergyDTA,
172 incidentParticle, targetParticle,
173 atomicWeight, atomicNumber);
174
175 if (!successful)
177 excitationEnergyGNP, excitationEnergyDTA,
178 incidentParticle, targetParticle,
179 atomicWeight, atomicNumber);
180 if (!successful)
181 QuasiElasticScattering(successful, pv, vecLength,
182 excitationEnergyGNP, excitationEnergyDTA,
183 incidentParticle, targetParticle,
184 atomicWeight, atomicNumber);
185 if (!successful)
186 ElasticScattering(successful, pv, vecLength,
187 incidentParticle,
188 atomicWeight, atomicNumber);
189
190 if (!successful)
191 G4cout << "GHEInelasticInteraction::ApplyYourself fails to produce final state particles"
192 << G4endl;
193
195 delete [] pv;
197 return &theParticleChange;
198}
@ stopAndKill
double G4double
Definition: G4Types.hh:64
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
#define G4UniformRand()
Definition: Randomize.hh:53
void MediumEnergyClusterProduction(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
void ElasticScattering(G4bool &successful, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, G4double atomicWeight, G4double atomicNumber)
void QuasiElasticScattering(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
void FillParticleChange(G4HEVector pv[], G4int aVecLength)
void HighEnergyClusterProduction(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
G4double NuclearExcitation(G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber, G4double &excitationEnergyCascade, G4double &excitationEnergyEvaporation)
void MediumEnergyCascading(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
G4double NuclearInelasticity(G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber)
void StrangeParticlePairProduction(const G4double availableEnergy, const G4double centerOfMassEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle)
void HighEnergyCascading(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
void FirstIntInCasSigmaMinus(G4bool &inElastic, const G4double availableEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, const G4double atomicWeight)
G4double getMass() const
Definition: G4HEVector.cc:361
void setDefinition(G4String name)
Definition: G4HEVector.cc:812
void SetStatusChange(G4HadFinalStateStatus aS)
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115

◆ FirstIntInCasSigmaMinus()

void G4HESigmaMinusInelastic::FirstIntInCasSigmaMinus ( G4bool inElastic,
const G4double  availableEnergy,
G4HEVector  pv[],
G4int vecLen,
const G4HEVector incidentParticle,
const G4HEVector targetParticle,
const G4double  atomicWeight 
)

Definition at line 202 of file G4HESigmaMinusInelastic.cc.

217{
218 static const G4double expxu = 82.; // upper bound for arg. of exp
219 static const G4double expxl = -expxu; // lower bound for arg. of exp
220
221 static const G4double protb = 0.7;
222 static const G4double neutb = 0.7;
223 static const G4double c = 1.25;
224
225 static const G4int numMul = 1200;
226 static const G4int numSec = 60;
227
229 G4int protonCode = Proton.getCode();
230
231 G4int targetCode = targetParticle.getCode();
232 G4double incidentTotalMomentum = incidentParticle.getTotalMomentum();
233
234 static G4bool first = true;
235 static G4double protmul[numMul], protnorm[numSec]; // proton constants
236 static G4double neutmul[numMul], neutnorm[numSec]; // neutron constants
237
238 // misc. local variables
239 // npos = number of pi+, nneg = number of pi-, nzero = number of pi0
240
241 G4int i, counter, nt, npos, nneg, nzero;
242
243 if (first) { // computation of normalization constants will only be done once
244 first = false;
245 for( i=0; i<numMul; i++ )protmul[i] = 0.0;
246 for( i=0; i<numSec; i++ )protnorm[i] = 0.0;
247 counter = -1;
248 for (npos = 0; npos < (numSec/3); npos++) {
249 for (nneg = std::max(0,npos-1); nneg <= (npos+1); nneg++) {
250 for (nzero = 0; nzero < numSec/3; nzero++) {
251 if (++counter < numMul) {
252 nt = npos+nneg+nzero;
253 if ((nt > 0) && (nt<=numSec) ) {
254 protmul[counter] = pmltpc(npos,nneg,nzero,nt,protb,c);
255 protnorm[nt-1] += protmul[counter];
256 }
257 }
258 }
259 }
260 }
261
262 for (i = 0; i < numMul; i++) neutmul[i] = 0.0;
263 for (i = 0; i < numSec; i++) neutnorm[i] = 0.0;
264 counter = -1;
265 for (npos = 0; npos < numSec/3; npos++) {
266 for (nneg = npos; nneg <= (npos+2); nneg++) {
267 for (nzero = 0; nzero < numSec/3; nzero++) {
268 if (++counter < numMul) {
269 nt = npos+nneg+nzero;
270 if ((nt>0) && (nt<=numSec) ) {
271 neutmul[counter] = pmltpc(npos,nneg,nzero,nt,neutb,c);
272 neutnorm[nt-1] += neutmul[counter];
273 }
274 }
275 }
276 }
277 }
278 for (i = 0; i < numSec; i++) {
279 if (protnorm[i] > 0.0 )protnorm[i] = 1.0/protnorm[i];
280 if (neutnorm[i] > 0.0 )neutnorm[i] = 1.0/neutnorm[i];
281 }
282 } // end of initialization
283
284 pv[0] = incidentParticle; // initialize the first two places
285 pv[1] = targetParticle; // the same as beam and target
286 vecLen = 2;
287
288 if (!inElastic) { // quasi-elastic scattering, no pions produced
289 G4double cech[] = {0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.06, 0.04, 0.005, 0.};
290 G4int iplab = G4int( std::min( 9.0, incidentTotalMomentum*2.5 ) );
291 if (G4UniformRand() < cech[iplab]/std::pow(atomicWeight,0.42) ) {
292 G4double ran = G4UniformRand();
293 if (targetCode == neutronCode) {
294 pv[0] = Neutron;
295 pv[1] = SigmaMinus;
296 } else {
297 if (ran < 0.2) {
298 pv[0] = SigmaZero;
299 pv[1] = Neutron;
300 } else if (ran < 0.4) {
301 pv[0] = Lambda;
302 pv[1] = Neutron;
303 } else if (ran < 0.6) {
304 pv[0] = Proton;
305 pv[1] = SigmaMinus;
306 } else if (ran < 0.8) {
307 pv[0] = Neutron;
308 pv[1] = SigmaZero;
309 } else {
310 pv[0] = Neutron;
311 pv[1] = Lambda;
312 }
313 }
314 }
315 return;
316
317 } else if (availableEnergy <= PionPlus.getMass()) return;
318
319 // inelastic scattering
320 npos = 0, nneg = 0, nzero = 0;
321
322 // number of total particles vs. centre of mass Energy - 2*proton mass
323 G4double aleab = std::log(availableEnergy);
324 G4double n = 3.62567+aleab*(0.665843+aleab*(0.336514
325 + aleab*(0.117712+0.0136912*aleab))) - 2.0;
326
327 // normalization constant for kno-distribution.
328 // calculate first the sum of all constants, check for numerical problems.
329 G4double test, dum, anpn = 0.0;
330
331 for (nt = 1; nt <= numSec; nt++) {
332 test = std::exp(std::min(expxu, std::max(expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
333 dum = pi*nt/(2.0*n*n);
334 if (std::fabs(dum) < 1.0) {
335 if (test >= 1.0e-10) anpn += dum*test;
336 } else {
337 anpn += dum*test;
338 }
339 }
340
341 G4double ran = G4UniformRand();
342 G4double excs = 0.0;
343 if (targetCode == protonCode) {
344 counter = -1;
345 for (npos=0; npos<numSec/3; npos++) {
346 for (nneg=std::max(0,npos-1); nneg<=(npos+1); nneg++) {
347 for (nzero=0; nzero<numSec/3; nzero++) {
348 if (++counter < numMul) {
349 nt = npos+nneg+nzero;
350 if ( (nt>0) && (nt<=numSec) ) {
351 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
352 dum = (pi/anpn)*nt*protmul[counter]*protnorm[nt-1]/(2.0*n*n);
353 if (std::fabs(dum) < 1.0) {
354 if( test >= 1.0e-10 )excs += dum*test;
355 } else {
356 excs += dum*test;
357 }
358 if (ran < excs) goto outOfLoop; //----------------------->
359 }
360 }
361 }
362 }
363 }
364
365 // 3 previous loops continued to the end
366 inElastic = false; // quasi-elastic scattering
367 return;
368
369 } else { // target must be a neutron
370 counter = -1;
371 for (npos=0; npos<numSec/3; npos++) {
372 for (nneg=npos; nneg<=(npos+2); nneg++) {
373 for (nzero=0; nzero<numSec/3; nzero++) {
374 if (++counter < numMul) {
375 nt = npos+nneg+nzero;
376 if ( (nt>=1) && (nt<=numSec) ) {
377 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
378 dum = (pi/anpn)*nt*neutmul[counter]*neutnorm[nt-1]/(2.0*n*n);
379 if (std::fabs(dum) < 1.0) {
380 if( test >= 1.0e-10 )excs += dum*test;
381 } else {
382 excs += dum*test;
383 }
384 if (ran < excs) goto outOfLoop; // -------------------->
385 }
386 }
387 }
388 }
389 }
390 // 3 previous loops continued to the end
391
392 inElastic = false; // quasi-elastic scattering.
393 return;
394 }
395
396 outOfLoop: // <-------------------------------------------
397
398 ran = G4UniformRand();
399 if (targetCode == neutronCode) {
400 if( npos == nneg)
401 {
402 }
403 else if (npos == (nneg-1))
404 {
405 if( ran < 0.25)
406 {
407 pv[0] = SigmaZero;
408 }
409 else if (ran < 0.5)
410 {
411 pv[0] = Lambda;
412 }
413 else
414 {
415 pv[1] = Proton;
416 }
417 }
418 else
419 {
420 if(ran < 0.5)
421 {
422 pv[0] = SigmaZero;
423 pv[1] = Proton;
424 }
425 else
426 {
427 pv[0] = Lambda;
428 pv[1] = Proton;
429 }
430 }
431 } else {
432 if (npos == nneg)
433 {
434 if (ran < 0.5)
435 {
436 }
437 else if (ran < 0.75)
438 {
439 pv[0] = SigmaZero;
440 pv[1] = Neutron;
441 }
442 else
443 {
444 pv[0] = Lambda;
445 pv[1] = Neutron;
446 }
447 }
448 else if (npos == (nneg+1))
449 {
450 pv[1] = Neutron;
451 }
452 else
453 {
454 if (ran < 0.5)
455 {
456 pv[0] = SigmaZero;
457 }
458 else
459 {
460 pv[0] = Lambda;
461 }
462 }
463 }
464
465 nt = npos + nneg + nzero;
466 while (nt > 0) {
467 G4double rnd = G4UniformRand();
468 if (rnd < (G4double)npos/nt) {
469 if (npos > 0) {
470 pv[vecLen++] = PionPlus;
471 npos--;
472 }
473 } else if (rnd < (G4double)(npos+nneg)/nt) {
474 if (nneg > 0) {
475 pv[vecLen++] = PionMinus;
476 nneg--;
477 }
478 } else {
479 if (nzero > 0) {
480 pv[vecLen++] = PionZero;
481 nzero--;
482 }
483 }
484 nt = npos + nneg + nzero;
485 }
486
487 if (verboseLevel > 1) {
488 G4cout << "Particles produced: " ;
489 G4cout << pv[0].getCode() << " " ;
490 G4cout << pv[1].getCode() << " " ;
491 for (i = 2; i < vecLen; i++) G4cout << pv[i].getCode() << " " ;
492 G4cout << G4endl;
493 }
494
495 return;
496}
@ neutronCode
G4HEVector PionPlus
G4HEVector SigmaZero
G4double pmltpc(G4int np, G4int nm, G4int nz, G4int n, G4double b, G4double c)
G4HEVector Lambda
G4HEVector SigmaMinus
G4HEVector Neutron
G4HEVector PionMinus
G4HEVector PionZero
G4HEVector Proton
G4int getCode() const
Definition: G4HEVector.cc:426
G4double getTotalMomentum() const
Definition: G4HEVector.cc:166
const G4double pi

Referenced by ApplyYourself().

◆ GetNumberOfSecondaries()

G4int G4HESigmaMinusInelastic::GetNumberOfSecondaries ( )
inline

Definition at line 65 of file G4HESigmaMinusInelastic.hh.

65{return vecLength;}

◆ ModelDescription()

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

Reimplemented from G4HadronicInteraction.

Definition at line 55 of file G4HESigmaMinusInelastic.cc.

56{
57 outFile << "G4HESigmaMinusInelastic is one of the High Energy Parameterized\n"
58 << "(HEP) models used to implement inelastic Sigma- scattering\n"
59 << "from nuclei. It is a re-engineered version of the GHEISHA\n"
60 << "code of H. Fesefeldt. It divides the initial collision\n"
61 << "products into backward- and forward-going clusters which are\n"
62 << "then decayed into final state hadrons. The model does not\n"
63 << "conserve energy on an event-by-event basis. It may be\n"
64 << "applied to sigmas with initial energies above 20 GeV.\n";
65}

Member Data Documentation

◆ vecLength

G4int G4HESigmaMinusInelastic::vecLength

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