Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
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