Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4HEAntiLambdaInelastic.cc
Go to the documentation of this file.
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// $Id$
27//
28
29// G4 Process: Gheisha High Energy Collision model.
30// This includes the high energy cascading model, the two-body-resonance model
31// and the low energy two-body model. Not included are the low energy stuff
32// like nuclear reactions, nuclear fission without any cascading and all
33// processes for particles at rest.
34// First work done by J.L.Chuma and F.W.Jones, TRIUMF, June 96.
35// H. Fesefeldt, RWTH-Aachen, 23-October-1996
36
37#include "G4HEAntiLambdaInelastic.hh"
38#include "globals.hh"
39#include "G4ios.hh"
40#include "G4PhysicalConstants.hh"
41#include "G4SystemOfUnits.hh"
42
43G4HEAntiLambdaInelastic::G4HEAntiLambdaInelastic(const G4String& name)
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 G4HEAntiLambdaInelastic is being deprecated and will\n"
52 << "disappear in Geant4 version 10.0" << G4endl;
53}
54
55
56void G4HEAntiLambdaInelastic::ModelDescription(std::ostream& outFile) const
57{
58 outFile << "G4HEAntiLambdaInelastic is one of the High Energy\n"
59 << "Parameterized (HEP) models used to implement inelastic\n"
60 << "anti-Lambda scattering from nuclei. It is a re-engineered\n"
61 << "version of the GHEISHA code of H. Fesefeldt. It divides the\n"
62 << "initial collision products into backward- and forward-going\n"
63 << "clusters which are then decayed into final state hadrons.\n"
64 << "The model does not conserve energy on an event-by-event\n"
65 << "basis. It may be applied to anti-Lambdas with initial energies\n"
66 << "above 20 GeV.\n";
67}
68
69
70G4HadFinalState*
71G4HEAntiLambdaInelastic::ApplyYourself(const G4HadProjectile &aTrack,
72 G4Nucleus &targetNucleus)
73{
74 G4HEVector* pv = new G4HEVector[MAXPART];
75 const G4HadProjectile *aParticle = &aTrack;
76 const G4double atomicWeight = targetNucleus.GetA_asInt();
77 const G4double atomicNumber = targetNucleus.GetZ_asInt();
78 G4HEVector incidentParticle(aParticle);
79
80 G4int incidentCode = incidentParticle.getCode();
81 G4double incidentMass = incidentParticle.getMass();
82 G4double incidentTotalEnergy = incidentParticle.getEnergy();
83
84 G4double incidentKineticEnergy = incidentTotalEnergy - incidentMass;
85
86 if (incidentKineticEnergy < 1.)
87 G4cout << "GHEAntiLambdaInelastic: incident energy < 1 GeV" << G4endl;
88
89 if (verboseLevel > 1) {
90 G4cout << "G4HEAntiLambdaInelastic::ApplyYourself" << G4endl;
91 G4cout << "incident particle " << incidentParticle.getName()
92 << "mass " << incidentMass
93 << "kinetic energy " << incidentKineticEnergy
94 << G4endl;
95 G4cout << "target material with (A,Z) = ("
96 << atomicWeight << "," << atomicNumber << ")" << G4endl;
97 }
98
99 G4double inelasticity = NuclearInelasticity(incidentKineticEnergy,
100 atomicWeight, atomicNumber);
101 if (verboseLevel > 1)
102 G4cout << "nuclear inelasticity = " << inelasticity << G4endl;
103
104 incidentKineticEnergy -= inelasticity;
105
106 G4double excitationEnergyGNP = 0.;
107 G4double excitationEnergyDTA = 0.;
108
109 G4double excitation = NuclearExcitation(incidentKineticEnergy,
110 atomicWeight, atomicNumber,
111 excitationEnergyGNP,
112 excitationEnergyDTA);
113 if (verboseLevel > 1)
114 G4cout << "nuclear excitation = " << excitation << excitationEnergyGNP
115 << excitationEnergyDTA << G4endl;
116
117 incidentKineticEnergy -= excitation;
118 incidentTotalEnergy = incidentKineticEnergy + incidentMass;
119 // incidentTotalMomentum = std::sqrt( (incidentTotalEnergy-incidentMass)
120 // *(incidentTotalEnergy+incidentMass));
121 // DHW 19 May 2011: variable set but not used
122
123 G4HEVector targetParticle;
124 if (G4UniformRand() < atomicNumber/atomicWeight) {
125 targetParticle.setDefinition("Proton");
126 } else {
127 targetParticle.setDefinition("Neutron");
128 }
129
130 G4double targetMass = targetParticle.getMass();
131 G4double centerOfMassEnergy =
132 std::sqrt( incidentMass*incidentMass + targetMass*targetMass
133 + 2.0*targetMass*incidentTotalEnergy);
134 G4double availableEnergy = centerOfMassEnergy - targetMass - incidentMass;
135
136 G4bool inElastic = true;
137 vecLength = 0;
138
139 if (verboseLevel > 1)
140 G4cout << "ApplyYourself: CallFirstIntInCascade for particle "
141 << incidentCode << G4endl;
142
143 G4bool successful = false;
144
145 FirstIntInCasAntiLambda(inElastic, availableEnergy, pv, vecLength,
146 incidentParticle, targetParticle, atomicWeight);
147
148 if (verboseLevel > 1)
149 G4cout << "ApplyYourself::StrangeParticlePairProduction" << G4endl;
150
151 if ((vecLength > 0) && (availableEnergy > 1.))
152 StrangeParticlePairProduction(availableEnergy, centerOfMassEnergy,
153 pv, vecLength,
154 incidentParticle, targetParticle);
155 HighEnergyCascading(successful, pv, vecLength,
156 excitationEnergyGNP, excitationEnergyDTA,
157 incidentParticle, targetParticle,
158 atomicWeight, atomicNumber);
159 if (!successful)
160 HighEnergyClusterProduction(successful, pv, vecLength,
161 excitationEnergyGNP, excitationEnergyDTA,
162 incidentParticle, targetParticle,
163 atomicWeight, atomicNumber);
164 if (!successful)
165 MediumEnergyCascading(successful, pv, vecLength,
166 excitationEnergyGNP, excitationEnergyDTA,
167 incidentParticle, targetParticle,
168 atomicWeight, atomicNumber);
169
170 if (!successful)
171 MediumEnergyClusterProduction(successful, pv, vecLength,
172 excitationEnergyGNP, excitationEnergyDTA,
173 incidentParticle, targetParticle,
174 atomicWeight, atomicNumber);
175 if (!successful)
176 QuasiElasticScattering(successful, pv, vecLength,
177 excitationEnergyGNP, excitationEnergyDTA,
178 incidentParticle, targetParticle,
179 atomicWeight, atomicNumber);
180 if (!successful)
181 ElasticScattering(successful, pv, vecLength,
182 incidentParticle,
183 atomicWeight, atomicNumber);
184
185 if (!successful)
186 G4cout << "GHEInelasticInteraction::ApplyYourself fails to produce final state particles"
187 << G4endl;
188
189 FillParticleChange(pv, vecLength);
190 delete [] pv;
191 theParticleChange.SetStatusChange(stopAndKill);
192 return & theParticleChange;
193}
194
195
196void
197G4HEAntiLambdaInelastic::FirstIntInCasAntiLambda(G4bool& inElastic,
198 const G4double availableEnergy,
199 G4HEVector pv[],
200 G4int &vecLen,
201 const G4HEVector& incidentParticle,
202 const G4HEVector& targetParticle,
203 const G4double atomicWeight)
204
205// AntiLambda undergoes interaction with nucleon within a nucleus.
206// Check if it is energetically possible to produce pions/kaons. If not,
207// assume nuclear excitation occurs and input particle is degraded in
208// energy. No other particles are produced.
209// If reaction is possible, find the correct number of pions/protons/neutrons
210// produced using an interpolation to multiplicity data. Replace some pions or
211// protons/neutrons by kaons or strange baryons according to the average
212// multiplicity per inelastic reaction.
213{
214 static const G4double expxu = 82.; // upper bound for arg. of exp
215 static const G4double expxl = -expxu; // lower bound for arg. of exp
216
217 static const G4double protb = 0.7;
218 static const G4double neutb = 0.7;
219 static const G4double c = 1.25;
220
221 static const G4int numMul = 1200;
222 static const G4int numMulAn = 400;
223 static const G4int numSec = 60;
224
225 G4int protonCode = Proton.getCode();
226
227 G4int targetCode = targetParticle.getCode();
228 G4double incidentTotalMomentum = incidentParticle.getTotalMomentum();
229
230 static G4bool first = true;
231 static G4double protmul[numMul], protnorm[numSec]; // proton constants
232 static G4double protmulAn[numMulAn],protnormAn[numSec];
233 static G4double neutmul[numMul], neutnorm[numSec]; // neutron constants
234 static G4double neutmulAn[numMulAn],neutnormAn[numSec];
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 )
242 { // compute normalization constants, this will only be done once
243 first = false;
244 for( i=0; i<numMul ; i++ ) protmul[i] = 0.0;
245 for( i=0; i<numSec ; i++ ) protnorm[i] = 0.0;
246 counter = -1;
247 for( npos=0; npos<(numSec/3); npos++ )
248 {
249 for( nneg=std::max(0,npos-2); nneg<=(npos+1); nneg++ )
250 {
251 for( nzero=0; nzero<numSec/3; nzero++ )
252 {
253 if( ++counter < numMul )
254 {
255 nt = npos+nneg+nzero;
256 if( (nt>0) && (nt<=numSec) )
257 {
258 protmul[counter] = pmltpc(npos,nneg,nzero,nt,protb,c);
259 protnorm[nt-1] += protmul[counter];
260 }
261 }
262 }
263 }
264 }
265 for( i=0; i<numMul; i++ )neutmul[i] = 0.0;
266 for( i=0; i<numSec; i++ )neutnorm[i] = 0.0;
267 counter = -1;
268 for( npos=0; npos<numSec/3; npos++ )
269 {
270 for( nneg=std::max(0,npos-1); nneg<=(npos+2); nneg++ )
271 {
272 for( nzero=0; nzero<numSec/3; nzero++ )
273 {
274 if( ++counter < numMul )
275 {
276 nt = npos+nneg+nzero;
277 if( (nt>0) && (nt<=numSec) )
278 {
279 neutmul[counter] = pmltpc(npos,nneg,nzero,nt,neutb,c);
280 neutnorm[nt-1] += neutmul[counter];
281 }
282 }
283 }
284 }
285 }
286 for( i=0; i<numSec; i++ )
287 {
288 if( protnorm[i] > 0.0 )protnorm[i] = 1.0/protnorm[i];
289 if( neutnorm[i] > 0.0 )neutnorm[i] = 1.0/neutnorm[i];
290 }
291 // annihilation
292 for( i=0; i<numMulAn ; i++ ) protmulAn[i] = 0.0;
293 for( i=0; i<numSec ; i++ ) protnormAn[i] = 0.0;
294 counter = -1;
295 for( npos=1; npos<(numSec/3); npos++ )
296 {
297 nneg = std::max(0,npos-1);
298 for( nzero=0; nzero<numSec/3; nzero++ )
299 {
300 if( ++counter < numMulAn )
301 {
302 nt = npos+nneg+nzero;
303 if( (nt>1) && (nt<=numSec) )
304 {
305 protmulAn[counter] = pmltpc(npos,nneg,nzero,nt,protb,c);
306 protnormAn[nt-1] += protmulAn[counter];
307 }
308 }
309 }
310 }
311 for( i=0; i<numMulAn; i++ ) neutmulAn[i] = 0.0;
312 for( i=0; i<numSec; i++ ) neutnormAn[i] = 0.0;
313 counter = -1;
314 for( npos=0; npos<numSec/3; npos++ )
315 {
316 nneg = npos;
317 for( nzero=0; nzero<numSec/3; nzero++ )
318 {
319 if( ++counter < numMulAn )
320 {
321 nt = npos+nneg+nzero;
322 if( (nt>1) && (nt<=numSec) )
323 {
324 neutmulAn[counter] = pmltpc(npos,nneg,nzero,nt,neutb,c);
325 neutnormAn[nt-1] += neutmulAn[counter];
326 }
327 }
328 }
329 }
330 for( i=0; i<numSec; i++ )
331 {
332 if( protnormAn[i] > 0.0 )protnormAn[i] = 1.0/protnormAn[i];
333 if( neutnormAn[i] > 0.0 )neutnormAn[i] = 1.0/neutnormAn[i];
334 }
335 } // end of initialization
336
337
338 // initialize the first two places
339 // the same as beam and target
340 pv[0] = incidentParticle;
341 pv[1] = targetParticle;
342 vecLen = 2;
343
344 if( !inElastic )
345 { // some two-body reactions
346 G4double cech[] = {0.50, 0.45, 0.40, 0.35, 0.30, 0.25, 0.06, 0.04, 0.005, 0.};
347
348 G4int iplab = std::min(9, G4int( incidentTotalMomentum*2.5));
349 if( G4UniformRand() < cech[iplab]/std::pow(atomicWeight,0.42) )
350 {
351 G4double ran = G4UniformRand();
352
353 if ( targetCode == protonCode)
354 {
355 if(ran < 0.2)
356 {
357 pv[0] = AntiSigmaZero;
358 }
359 else if (ran < 0.4)
360 {
361 pv[0] = AntiSigmaMinus;
362 pv[1] = Neutron;
363 }
364 else if (ran < 0.6)
365 {
366 pv[0] = Proton;
367 pv[1] = AntiLambda;
368 }
369 else if (ran < 0.8)
370 {
371 pv[0] = Proton;
372 pv[1] = AntiSigmaZero;
373 }
374 else
375 {
376 pv[0] = Neutron;
377 pv[1] = AntiSigmaMinus;
378 }
379 }
380 else
381 {
382 if (ran < 0.2)
383 {
384 pv[0] = AntiSigmaZero;
385 }
386 else if (ran < 0.4)
387 {
388 pv[0] = AntiSigmaPlus;
389 pv[1] = Proton;
390 }
391 else if (ran < 0.6)
392 {
393 pv[0] = Neutron;
394 pv[1] = AntiLambda;
395 }
396 else if (ran < 0.8)
397 {
398 pv[0] = Neutron;
399 pv[1] = AntiSigmaZero;
400 }
401 else
402 {
403 pv[0] = Proton;
404 pv[1] = AntiSigmaPlus;
405 }
406 }
407 }
408 return;
409 }
410 else if (availableEnergy <= PionPlus.getMass())
411 return;
412
413 // inelastic scattering
414
415 npos = 0; nneg = 0; nzero = 0;
416 G4double anhl[] = {1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 1.00, 0.97, 0.88,
417 0.85, 0.81, 0.75, 0.64, 0.64, 0.55, 0.55, 0.45, 0.47, 0.40,
418 0.39, 0.36, 0.33, 0.10, 0.01};
419 G4int iplab = G4int( incidentTotalMomentum*10.);
420 if ( iplab > 9) iplab = 10 + G4int( (incidentTotalMomentum -1.)*5. );
421 if ( iplab > 14) iplab = 15 + G4int( incidentTotalMomentum -2. );
422 if ( iplab > 22) iplab = 23 + G4int( (incidentTotalMomentum -10.)/10.);
423 iplab = std::min(24, iplab);
424
425 if (G4UniformRand() > anhl[iplab]) { // non- annihilation channels
426
427 // number of total particles vs. centre of mass Energy - 2*proton mass
428 G4double aleab = std::log(availableEnergy);
429 G4double n = 3.62567+aleab*(0.665843+aleab*(0.336514
430 + aleab*(0.117712+0.0136912*aleab))) - 2.0;
431
432 // normalization constant for kno-distribution.
433 // calculate first the sum of all constants, check for numerical problems.
434 G4double test, dum, anpn = 0.0;
435
436 for (nt = 1; nt <= numSec; nt++) {
437 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
438 dum = pi*nt/(2.0*n*n);
439 if (std::fabs(dum) < 1.0) {
440 if (test >= 1.0e-10) anpn += dum*test;
441 } else {
442 anpn += dum*test;
443 }
444 }
445
446 G4double ran = G4UniformRand();
447 G4double excs = 0.0;
448 if (targetCode == protonCode) {
449 counter = -1;
450 for (npos = 0; npos < numSec/3; npos++) {
451 for (nneg = std::max(0,npos-2); nneg <= (npos+1); nneg++) {
452 for (nzero = 0; nzero < numSec/3; nzero++) {
453 if (++counter < numMul) {
454 nt = npos+nneg+nzero;
455 if ((nt > 0) && (nt <= numSec) ) {
456 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
457 dum = (pi/anpn)*nt*protmul[counter]*protnorm[nt-1]/(2.0*n*n);
458
459 if (std::fabs(dum) < 1.0) {
460 if (test >= 1.0e-10) excs += dum*test;
461 } else {
462 excs += dum*test;
463 }
464
465 if (ran < excs) goto outOfLoop; //----------------------->
466 }
467 }
468 }
469 }
470 }
471 // 3 previous loops continued to the end
472 inElastic = false; // quasi-elastic scattering
473 return;
474 } else { // target must be a neutron
475 counter = -1;
476 for( npos=0; npos<numSec/3; npos++ )
477 {
478 for( nneg=std::max(0,npos-1); nneg<=(npos+2); nneg++ )
479 {
480 for( nzero=0; nzero<numSec/3; nzero++ )
481 {
482 if( ++counter < numMul )
483 {
484 nt = npos+nneg+nzero;
485 if( (nt>0) && (nt<=numSec) )
486 {
487 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
488 dum = (pi/anpn)*nt*neutmul[counter]*neutnorm[nt-1]/(2.0*n*n);
489 if (std::fabs(dum) < 1.0) {
490 if( test >= 1.0e-10 )excs += dum*test;
491 } else {
492 excs += dum*test;
493 }
494
495 if (ran < excs) goto outOfLoop; // -------------------------->
496 }
497 }
498 }
499 }
500 }
501 // 3 previous loops continued to the end
502 inElastic = false; // quasi-elastic scattering.
503 return;
504 }
505
506 outOfLoop: // <------------------------------------------------------------------------
507
508 ran = G4UniformRand();
509
510 if( targetCode == protonCode)
511 {
512 if( npos == nneg)
513 {
514 if (ran < 0.40)
515 {
516 }
517 else if (ran < 0.8)
518 {
519 pv[0] = AntiSigmaZero;
520 }
521 else
522 {
523 pv[0] = AntiSigmaMinus;
524 pv[1] = Neutron;
525 }
526 }
527 else if (npos == (nneg+1))
528 {
529 if( ran < 0.25)
530 {
531 pv[1] = Neutron;
532 }
533 else if (ran < 0.5)
534 {
535 pv[0] = AntiSigmaZero;
536 pv[1] = Neutron;
537 }
538 else
539 {
540 pv[0] = AntiSigmaPlus;
541 }
542 }
543 else if (npos == (nneg-1))
544 {
545 pv[0] = AntiSigmaMinus;
546 }
547 else
548 {
549 pv[0] = AntiSigmaPlus;
550 pv[1] = Neutron;
551 }
552 }
553 else
554 {
555 if( npos == nneg)
556 {
557 if (ran < 0.4)
558 {
559 }
560 else if(ran < 0.8)
561 {
562 pv[0] = AntiSigmaZero;
563 }
564 else
565 {
566 pv[0] = AntiSigmaPlus;
567 pv[1] = Proton;
568 }
569 }
570 else if ( npos == (nneg-1))
571 {
572 if (ran < 0.5)
573 {
574 pv[0] = AntiSigmaMinus;
575 }
576 else if (ran < 0.75)
577 {
578 pv[1] = Proton;
579 }
580 else
581 {
582 pv[0] = AntiSigmaZero;
583 pv[1] = Proton;
584 }
585 }
586 else if (npos == (nneg+1))
587 {
588 pv[0] = AntiSigmaPlus;
589 }
590 else
591 {
592 pv[0] = AntiSigmaMinus;
593 pv[1] = Proton;
594 }
595 }
596
597 }
598 else // annihilation
599 {
600 if ( availableEnergy > 2. * PionPlus.getMass() )
601 {
602
603 G4double aleab = std::log(availableEnergy);
604 G4double n = 3.62567+aleab*(0.665843+aleab*(0.336514
605 + aleab*(0.117712+0.0136912*aleab))) - 2.0;
606
607 // normalization constant for kno-distribution.
608 // calculate first the sum of all constants, check for numerical problems.
609 G4double test, dum, anpn = 0.0;
610
611 for (nt=2; nt<=numSec; nt++) {
612 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
613 dum = pi*nt/(2.0*n*n);
614
615 if (std::fabs(dum) < 1.0) {
616 if( test >= 1.0e-10 )anpn += dum*test;
617 } else {
618 anpn += dum*test;
619 }
620 }
621
622 G4double ran = G4UniformRand();
623 G4double excs = 0.0;
624 if (targetCode == protonCode) {
625 counter = -1;
626 for (npos=1; npos<numSec/3; npos++) {
627 nneg = npos-1;
628 for( nzero=0; nzero<numSec/3; nzero++ )
629 {
630 if( ++counter < numMulAn )
631 {
632 nt = npos+nneg+nzero;
633 if( (nt>1) && (nt<=numSec) ) {
634 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
635 dum = (pi/anpn)*nt*protmulAn[counter]*protnormAn[nt-1]/(2.0*n*n);
636
637 if (std::fabs(dum) < 1.0) {
638 if( test >= 1.0e-10 )excs += dum*test;
639 } else {
640 excs += dum*test;
641 }
642
643 if (ran < excs) goto outOfLoopAn; //----------------------->
644 }
645 }
646 }
647 }
648 // 3 previous loops continued to the end
649 inElastic = false; // quasi-elastic scattering
650 return;
651
652 } else { // target must be a neutron
653 counter = -1;
654 for (npos=0; npos<numSec/3; npos++) {
655 nneg = npos;
656 for( nzero=0; nzero<numSec/3; nzero++ ) {
657 if (++counter < numMulAn) {
658 nt = npos+nneg+nzero;
659 if( (nt>1) && (nt<=numSec) ) {
660 test = std::exp( std::min( expxu, std::max( expxl, -(pi/4.0)*(nt*nt)/(n*n) ) ) );
661 dum = (pi/anpn)*nt*neutmulAn[counter]*neutnormAn[nt-1]/(2.0*n*n);
662
663 if (std::fabs(dum) < 1.0) {
664 if( test >= 1.0e-10 )excs += dum*test;
665 } else {
666 excs += dum*test;
667 }
668
669 if (ran < excs) goto outOfLoopAn; // -------------------------->
670 }
671 }
672 }
673 }
674
675 inElastic = false; // quasi-elastic scattering.
676 return;
677 }
678 outOfLoopAn: // <---------------------------------------------------------
679 vecLen = 0;
680 }
681 }
682
683 nt = npos + nneg + nzero;
684 while ( nt > 0)
685 {
686 G4double ran = G4UniformRand();
687 if ( ran < (G4double)npos/nt)
688 {
689 if( npos > 0 )
690 { pv[vecLen++] = PionPlus;
691 npos--;
692 }
693 }
694 else if ( ran < (G4double)(npos+nneg)/nt)
695 {
696 if( nneg > 0 )
697 {
698 pv[vecLen++] = PionMinus;
699 nneg--;
700 }
701 }
702 else
703 {
704 if( nzero > 0 )
705 {
706 pv[vecLen++] = PionZero;
707 nzero--;
708 }
709 }
710 nt = npos + nneg + nzero;
711 }
712 if (verboseLevel > 1)
713 {
714 G4cout << "Particles produced: " ;
715 G4cout << pv[0].getName() << " " ;
716 G4cout << pv[1].getName() << " " ;
717 for (i=2; i < vecLen; i++)
718 {
719 G4cout << pv[i].getName() << " " ;
720 }
721 G4cout << G4endl;
722 }
723 return;
724 }
725
726
727
728
729
730
731
732
733
734
stopAndKill
@ stopAndKill
Definition: G4HadFinalState.hh:42
G4double
double G4double
Definition: G4Types.hh:64
G4int
int G4int
Definition: G4Types.hh:66
G4bool
bool G4bool
Definition: G4Types.hh:67
G4endl
#define G4endl
Definition: G4ios.hh:52
G4cout
G4DLLIMPORT std::ostream G4cout
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:53
G4HEAntiLambdaInelastic::ModelDescription
virtual void ModelDescription(std::ostream &) const
Definition: G4HEAntiLambdaInelastic.cc:56
G4HEAntiLambdaInelastic::G4HEAntiLambdaInelastic
G4HEAntiLambdaInelastic(const G4String &name="G4HEAntiLambdaInelastic")
Definition: G4HEAntiLambdaInelastic.cc:43
G4HEAntiLambdaInelastic::FirstIntInCasAntiLambda
void FirstIntInCasAntiLambda(G4bool &inElastic, const G4double availableEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, const G4double atomicWeight)
Definition: G4HEAntiLambdaInelastic.cc:197
G4HEAntiLambdaInelastic::ApplyYourself
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
Definition: G4HEAntiLambdaInelastic.cc:71
G4HEInelastic::PionPlus
G4HEVector PionPlus
Definition: G4HEInelastic.hh:215
G4HEInelastic::AntiSigmaZero
G4HEVector AntiSigmaZero
Definition: G4HEInelastic.hh:234
G4HEInelastic::pmltpc
G4double pmltpc(G4int np, G4int nm, G4int nz, G4int n, G4double b, G4double c)
Definition: G4HEInelastic.cc:233
G4HEInelastic::AntiSigmaPlus
G4HEVector AntiSigmaPlus
Definition: G4HEInelastic.hh:233
G4HEInelastic::MediumEnergyClusterProduction
void MediumEnergyClusterProduction(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:4036
G4HEInelastic::ElasticScattering
void ElasticScattering(G4bool &successful, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:5183
G4HEInelastic::QuasiElasticScattering
void QuasiElasticScattering(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:4892
G4HEInelastic::Neutron
G4HEVector Neutron
Definition: G4HEInelastic.hh:226
G4HEInelastic::FillParticleChange
void FillParticleChange(G4HEVector pv[], G4int aVecLength)
Definition: G4HEInelastic.cc:57
G4HEInelastic::PionMinus
G4HEVector PionMinus
Definition: G4HEInelastic.hh:217
G4HEInelastic::HighEnergyClusterProduction
void HighEnergyClusterProduction(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:2097
G4HEInelastic::PionZero
G4HEVector PionZero
Definition: G4HEInelastic.hh:216
G4HEInelastic::NuclearExcitation
G4double NuclearExcitation(G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber, G4double &excitationEnergyCascade, G4double &excitationEnergyEvaporation)
Definition: G4HEInelastic.cc:179
G4HEInelastic::AntiSigmaMinus
G4HEVector AntiSigmaMinus
Definition: G4HEInelastic.hh:235
G4HEInelastic::AntiLambda
G4HEVector AntiLambda
Definition: G4HEInelastic.hh:229
G4HEInelastic::Proton
G4HEVector Proton
Definition: G4HEInelastic.hh:224
G4HEInelastic::MediumEnergyCascading
void MediumEnergyCascading(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:2970
G4HEInelastic::NuclearInelasticity
G4double NuclearInelasticity(G4double incidentKineticEnergy, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:149
G4HEInelastic::StrangeParticlePairProduction
void StrangeParticlePairProduction(const G4double availableEnergy, const G4double centerOfMassEnergy, G4HEVector pv[], G4int &vecLen, const G4HEVector &incidentParticle, const G4HEVector &targetParticle)
Definition: G4HEInelastic.cc:329
G4HEInelastic::HighEnergyCascading
void HighEnergyCascading(G4bool &successful, G4HEVector pv[], G4int &vecLen, G4double &excitationEnergyGNP, G4double &excitationEnergyDTA, const G4HEVector &incidentParticle, const G4HEVector &targetParticle, G4double atomicWeight, G4double atomicNumber)
Definition: G4HEInelastic.cc:598
G4HEVector::getEnergy
G4double getEnergy() const
Definition: G4HEVector.cc:313
G4HEVector::getMass
G4double getMass() const
Definition: G4HEVector.cc:361
G4HEVector::getCode
G4int getCode() const
Definition: G4HEVector.cc:426
G4HEVector::getTotalMomentum
G4double getTotalMomentum() const
Definition: G4HEVector.cc:166
G4HEVector::getName
G4String getName() const
Definition: G4HEVector.cc:431
G4HEVector::setDefinition
void setDefinition(G4String name)
Definition: G4HEVector.cc:812
G4HadFinalState::SetStatusChange
void SetStatusChange(G4HadFinalStateStatus aS)
Definition: G4HadFinalState.cc:81
G4Nucleus::GetA_asInt
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4Nucleus::GetZ_asInt
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115