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

#include <G4HELambdaInelastic.hh>

+ Inheritance diagram for G4HELambdaInelastic:

Public Member Functions

 G4HELambdaInelastic (const G4String &name="G4HELambdaInelastic")
 
 ~G4HELambdaInelastic ()
 
virtual void ModelDescription (std::ostream &) const
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4int GetNumberOfSecondaries ()
 
void FirstIntInCasLambda (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 52 of file G4HELambdaInelastic.hh.

Constructor & Destructor Documentation

◆ G4HELambdaInelastic()

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

Definition at line 43 of file G4HELambdaInelastic.cc.

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

◆ ~G4HELambdaInelastic()

G4HELambdaInelastic::~G4HELambdaInelastic ( )
inline

Definition at line 57 of file G4HELambdaInelastic.hh.

57{};

Member Function Documentation

◆ ApplyYourself()

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

Implements G4HadronicInteraction.

Definition at line 70 of file G4HELambdaInelastic.cc.

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

Referenced by G4HESigmaZeroInelastic::ApplyYourself().

◆ FirstIntInCasLambda()

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

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

Referenced by ApplyYourself().

◆ GetNumberOfSecondaries()

G4int G4HELambdaInelastic::GetNumberOfSecondaries ( )
inline

Definition at line 66 of file G4HELambdaInelastic.hh.

66{return vecLength;}

Referenced by G4HESigmaZeroInelastic::ApplyYourself().

◆ ModelDescription()

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

Reimplemented from G4HadronicInteraction.

Definition at line 56 of file G4HELambdaInelastic.cc.

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

Member Data Documentation

◆ vecLength

G4int G4HELambdaInelastic::vecLength

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