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

#include <G4HEOmegaMinusInelastic.hh>

+ Inheritance diagram for G4HEOmegaMinusInelastic:

Public Member Functions

 G4HEOmegaMinusInelastic ()
 
 ~G4HEOmegaMinusInelastic ()
 
virtual void ModelDescription (std::ostream &) const
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4int GetNumberOfSecondaries ()
 
void FirstIntInCasOmegaMinus (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 G4HEOmegaMinusInelastic.hh.

Constructor & Destructor Documentation

◆ G4HEOmegaMinusInelastic()

G4HEOmegaMinusInelastic::G4HEOmegaMinusInelastic ( )
inline

Definition at line 55 of file G4HEOmegaMinusInelastic.hh.

55 : G4HEInelastic("G4HEOmegaMinusInelastic")
56 {
57 vecLength = 0;
58 theMinEnergy = 20*CLHEP::GeV;
59 theMaxEnergy = 10*CLHEP::TeV;
60 MAXPART = 2048;
61 verboseLevel = 0;
62 G4cout << "WARNING: model G4HEOmegaMinusInelastic is being deprecated and will\n"
63 << "disappear in Geant4 version 10.0" << G4endl;
64 }
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout

◆ ~G4HEOmegaMinusInelastic()

G4HEOmegaMinusInelastic::~G4HEOmegaMinusInelastic ( )
inline

Definition at line 66 of file G4HEOmegaMinusInelastic.hh.

66{ };

Member Function Documentation

◆ ApplyYourself()

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

Implements G4HadronicInteraction.

Definition at line 57 of file G4HEOmegaMinusInelastic.cc.

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

◆ FirstIntInCasOmegaMinus()

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

Definition at line 190 of file G4HEOmegaMinusInelastic.cc.

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

Definition at line 75 of file G4HEOmegaMinusInelastic.hh.

75{return vecLength;}

◆ ModelDescription()

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

Reimplemented from G4HadronicInteraction.

Definition at line 42 of file G4HEOmegaMinusInelastic.cc.

43{
44 outFile << "G4HEOmegaMinusInelastic is one of the High Energy\n"
45 << "Parameterized (HEP) models used to implement inelastic\n"
46 << "Omega- scattering from nuclei. It is a re-engineered\n"
47 << "version of the GHEISHA code of H. Fesefeldt. It divides the\n"
48 << "initial collision products into backward- and forward-going\n"
49 << "clusters which are then decayed into final state hadrons.\n"
50 << "The model does not conserve energy on an event-by-event\n"
51 << "basis. It may be applied to Omega- with initial energies\n"
52 << "above 20 GeV.\n";
53}

Member Data Documentation

◆ vecLength

G4int G4HEOmegaMinusInelastic::vecLength

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