45 G4bool& incidentHasChanged,
66 if( vecLen == 0 )
return true;
70 if( currentParticle.
GetMass() == 0.0 || targetParticle.
GetMass() == 0.0 )
return true;
75 G4double centerofmassEnergy = std::sqrt( mOriginal*mOriginal +
76 targetMass*targetMass +
77 2.0*targetMass*etOriginal );
79 G4double availableEnergy = centerofmassEnergy-(targetMass+currentMass);
80 if( availableEnergy <= 1.0 )
return true;
104 const G4double avrs[] = {3.,4.,5.,6.,7.,8.,9.,10.,20.,30.,40.,50.};
112 ed <<
" While count exceeded " <<
G4endl;
113 while ((i<12) && (centerofmassEnergy>avrs[i]) ) {
137 i4 = i3 =
G4int(vecLen*rnd);
142 i4 =
G4int(vecLen*rnd);
148 const G4double avkkb[] = { 0.0015, 0.005, 0.012, 0.0285, 0.0525, 0.075,
149 0.0975, 0.123, 0.28, 0.398, 0.495, 0.573 };
150 const G4double avky[] = { 0.005, 0.03, 0.064, 0.095, 0.115, 0.13,
151 0.145, 0.155, 0.20, 0.205, 0.210, 0.212 };
152 const G4double avnnb[] = { 0.00001, 0.0001, 0.0006, 0.0025, 0.01, 0.02,
153 0.04, 0.05, 0.12, 0.15, 0.18, 0.20 };
155 avk = (
G4Log(avkkb[ibin])-
G4Log(avkkb[ibin-1]))*(centerofmassEnergy-avrs[ibin-1])
156 /(avrs[ibin]-avrs[ibin-1]) +
G4Log(avkkb[ibin-1]);
159 avy = (
G4Log(avky[ibin])-
G4Log(avky[ibin-1]))*(centerofmassEnergy-avrs[ibin-1])
160 /(avrs[ibin]-avrs[ibin-1]) +
G4Log(avky[ibin-1]);
163 avn = (
G4Log(avnnb[ibin])-
G4Log(avnnb[ibin-1]))*(centerofmassEnergy-avrs[ibin-1])
164 /(avrs[ibin]-avrs[ibin-1]) +
G4Log(avnnb[ibin-1]);
167 if( avk+avy+avn <= 0.0 )
return true;
169 if( currentMass < protonMass )avy /= 2.0;
170 if( targetMass < protonMass )avy = 0.0;
176 if( availableEnergy < 2.0 )
return true;
182 vec[0]->SetDefinition( aNeutron );
185 vec[0]->SetMayBeKilled(
false);
190 vec[0]->SetDefinition( aProton );
193 vec[0]->SetMayBeKilled(
false);
203 vec[i3]->SetDefinition( aNeutron );
204 vec[i4]->SetDefinition( anAntiNeutron );
205 vec[i3]->SetMayBeKilled(
false);
206 vec[i4]->SetMayBeKilled(
false);
210 vec[i3]->SetDefinition( aProton );
211 vec[i4]->SetDefinition( anAntiProton );
212 vec[i3]->SetMayBeKilled(
false);
213 vec[i4]->SetMayBeKilled(
false);
219 if( availableEnergy < 1.0 )
return true;
221 const G4double kkb[] = { 0.2500, 0.3750, 0.5000, 0.5625, 0.6250,
222 0.6875, 0.7500, 0.8750, 1.000 };
223 const G4int ipakkb1[] = { 10, 10, 10, 11, 11, 12, 12, 11, 12 };
224 const G4int ipakkb2[] = { 13, 11, 12, 11, 12, 11, 12, 13, 13 };
230 eda <<
" While count exceeded " <<
G4endl;
231 while( (i<9) && (ran>=kkb[i]) ) {
240 if( i == 9 )
return true;
248 vec[i3]->SetDefinition( aKaonPlus );
249 vec[i3]->SetMayBeKilled(
false);
252 vec[i3]->SetDefinition( aKaonZS );
253 vec[i3]->SetMayBeKilled(
false);
256 vec[i3]->SetDefinition( aKaonZL );
257 vec[i3]->SetMayBeKilled(
false);
288 vec[i4]->SetDefinition( aKaonZS );
289 vec[i4]->SetMayBeKilled(
false);
292 vec[i4]->SetDefinition( aKaonZL );
293 vec[i4]->SetMayBeKilled(
false);
296 vec[i4]->SetDefinition( aKaonMinus );
297 vec[i4]->SetMayBeKilled(
false);
302 }
else if( ran < avy ) {
303 if( availableEnergy < 1.6 )
return true;
305 const G4double ky[] = { 0.200, 0.300, 0.400, 0.550, 0.625, 0.700,
306 0.800, 0.850, 0.900, 0.950, 0.975, 1.000 };
307 const G4int ipaky1[] = { 18, 18, 18, 20, 20, 20, 21, 21, 21, 22, 22, 22 };
308 const G4int ipaky2[] = { 10, 11, 12, 10, 11, 12, 10, 11, 12, 10, 11, 12 };
309 const G4int ipakyb1[] = { 19, 19, 19, 23, 23, 23, 24, 24, 24, 25, 25, 25 };
310 const G4int ipakyb2[] = { 13, 12, 11, 13, 12, 11, 13, 12, 11, 13, 12, 11 };
316 edb <<
" While count exceeded " <<
G4endl;
317 while( (i<12) && (ran>ky[i]) ) {
326 if( i == 12 )
return true;
348 targetHasChanged =
true;
352 vec[i3]->SetDefinition( aKaonPlus );
353 vec[i3]->SetMayBeKilled(
false);
356 vec[i3]->SetDefinition( aKaonZS );
357 vec[i3]->SetMayBeKilled(
false);
360 vec[i3]->SetDefinition( aKaonZL );
361 vec[i3]->SetMayBeKilled(
false);
371 (currentMass > sigmaMinusMass) ) {
387 incidentHasChanged =
true;
391 vec[i3]->SetDefinition( aKaonZS );
392 vec[i3]->SetMayBeKilled(
false);
395 vec[i3]->SetDefinition( aKaonZL );
396 vec[i3]->SetMayBeKilled(
false);
399 vec[i3]->SetDefinition( aKaonMinus );
400 vec[i3]->SetMayBeKilled(
false);
420 incidentHasChanged =
true;
424 vec[i3]->SetDefinition( aKaonPlus );
425 vec[i3]->SetMayBeKilled(
false);
428 vec[i3]->SetDefinition( aKaonZS );
429 vec[i3]->SetMayBeKilled(
false);
432 vec[i3]->SetDefinition( aKaonZL );
433 vec[i3]->SetMayBeKilled(
false);
453 G4double energyCheck = centerofmassEnergy-(currentMass+targetMass);
454 for( i=0; i<vecLen; ++i )
456 energyCheck -= vec[i]->GetMass()/
GeV;
457 if( energyCheck < 0.0 )
461 for(j=i; j<vecLen; j++)
delete vec[j];
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
static G4Lambda * Lambda()
T max(const T t1, const T t2)
brief Return the largest of the two arguments
std::ostringstream G4ExceptionDescription
G4bool ReactionStage(const G4HadProjectile *, G4ReactionProduct &, G4bool &, const G4DynamicParticle *, G4ReactionProduct &, G4bool &, const G4Nucleus &, G4ReactionProduct &, G4FastVector< G4ReactionProduct, 256 > &, G4int &, G4bool, G4ReactionProduct &)
static G4AntiNeutron * AntiNeutron()
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
static G4AntiProton * AntiProton()
void SetDefinitionAndUpdateE(const G4ParticleDefinition *aParticleDefinition)
const G4ParticleDefinition * GetDefinition() const
static G4SigmaZero * SigmaZero()
static G4Proton * Proton()
static G4KaonMinus * KaonMinus()
static G4AntiSigmaMinus * AntiSigmaMinus()
G4double GetPDGMass() const
static G4KaonPlus * KaonPlus()
static G4SigmaMinus * SigmaMinus()
G4double G4Log(G4double x)
static G4SigmaPlus * SigmaPlus()
G4double GetTotalEnergy() const
static G4KaonZeroLong * KaonZeroLong()
G4ParticleDefinition * GetDefinition() const
void SetSide(const G4int sid)
static G4AntiSigmaZero * AntiSigmaZero()
void SetElement(G4int anIndex, Type *anElement)
static G4KaonZeroShort * KaonZeroShort()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
static G4Neutron * Neutron()
static G4AntiSigmaPlus * AntiSigmaPlus()
static constexpr double GeV
void SetMayBeKilled(const G4bool f)
static G4AntiLambda * AntiLambda()