Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | List of all members
G4INCL::BinaryCollisionAvatar Class Reference

#include <G4INCLBinaryCollisionAvatar.hh>

+ Inheritance diagram for G4INCL::BinaryCollisionAvatar:

Public Member Functions

 BinaryCollisionAvatar (G4double, G4double, G4INCL::Nucleus *, G4INCL::Particle *, G4INCL::Particle *)
 
virtual ~BinaryCollisionAvatar ()
 
G4INCL::IChannel * getChannel ()
 
ParticleList getParticles () const
 
virtual void preInteraction ()
 
virtual void postInteraction (FinalState *)
 
std::string dump () const
 
- Public Member Functions inherited from G4INCL::InteractionAvatar
 InteractionAvatar (G4double, G4INCL::Nucleus *, G4INCL::Particle *)
 
 InteractionAvatar (G4double, G4INCL::Nucleus *, G4INCL::Particle *, G4INCL::Particle *)
 
virtual ~InteractionAvatar ()
 
- Public Member Functions inherited from G4INCL::IAvatar
 IAvatar ()
 
 IAvatar (G4double time)
 
virtual ~IAvatar ()
 
FinalState * getFinalState ()
 
void fillFinalState (FinalState *fs)
 
G4double getTime () const
 
AvatarType getType () const
 
G4bool isACollision () const
 
G4bool isADecay () const
 
void setType (AvatarType t)
 
long getID () const
 
std::string toString ()
 

Static Public Member Functions

static void setCutNN (const G4double c)
 
static G4double getCutNN ()
 
static G4double getCutNNSquared ()
 
static G4double getBias ()
 Get the global bias factor.
 
static void setBias (const G4double b)
 Set the global bias factor.
 
- Static Public Member Functions inherited from G4INCL::InteractionAvatar
static void deleteBackupParticles ()
 Release the memory allocated for the backup particles.
 

Additional Inherited Members

- Static Public Attributes inherited from G4INCL::InteractionAvatar
static const G4double locEAccuracy = 1.E-4
 Target accuracy in the determination of the local-energy Q-value.
 
static const G4int maxIterLocE = 50
 Max number of iterations for the determination of the local-energy Q-value.
 
- Protected Member Functions inherited from G4INCL::InteractionAvatar
G4bool bringParticleInside (Particle *const p)
 
void preInteractionLocalEnergy (Particle *const p)
 Apply local-energy transformation, if appropriate.
 
void preInteractionBlocking ()
 Store the state of the particles before the interaction.
 
void preInteraction ()
 
void postInteraction (FinalState *)
 
void restoreParticles () const
 Restore the state of both particles.
 
G4bool shouldUseLocalEnergy () const
 true if the given avatar should use local energy
 
G4bool enforceEnergyConservation (FinalState *const fs)
 Enforce energy conservation.
 
- Protected Attributes inherited from G4INCL::InteractionAvatar
Nucleus * theNucleus
 
Particle * particle1
 
Particle * particle2
 
ThreeVector boostVector
 
G4double oldTotalEnergy
 
G4double oldXSec
 
G4bool isPiN
 
G4double weight
 
ParticleList modified
 
ParticleList created
 
ParticleList modifiedAndCreated
 
ParticleList Destroyed
 
ParticleList ModifiedAndDestroyed
 
- Protected Attributes inherited from G4INCL::IAvatar
G4double theTime
 
- Static Protected Attributes inherited from G4INCL::InteractionAvatar
static G4ThreadLocal Particle * backupParticle1 = NULL
 
static G4ThreadLocal Particle * backupParticle2 = NULL
 

Detailed Description

Definition at line 57 of file G4INCLBinaryCollisionAvatar.hh.

Constructor & Destructor Documentation

◆ BinaryCollisionAvatar()

G4INCL::BinaryCollisionAvatar::BinaryCollisionAvatar ( G4double time,
G4double crossSection,
G4INCL::Nucleus * n,
G4INCL::Particle * p1,
G4INCL::Particle * p2 )

Definition at line 131 of file G4INCLBinaryCollisionAvatar.cc.

133 : InteractionAvatar(time, n, p1, p2), theCrossSection(crossSection),
134 isParticle1Spectator(false),
135 isParticle2Spectator(false),
136 isElastic(false),
137 isStrangeProduction(false)
138 {
139 setType(CollisionAvatarType);
140 }
G4INCL::IAvatar::setType
void setType(AvatarType t)
Definition G4INCLIAvatar.hh:93
G4INCL::InteractionAvatar::InteractionAvatar
InteractionAvatar(G4double, G4INCL::Nucleus *, G4INCL::Particle *)
Definition G4INCLInteractionAvatar.cc:66
G4INCL::CollisionAvatarType
@ CollisionAvatarType
Definition G4INCLIAvatar.hh:67

◆ ~BinaryCollisionAvatar()

G4INCL::BinaryCollisionAvatar::~BinaryCollisionAvatar ( )
virtual

Definition at line 142 of file G4INCLBinaryCollisionAvatar.cc.

142 {
143 }

Member Function Documentation

◆ dump()

std::string G4INCL::BinaryCollisionAvatar::dump ( ) const
virtual

Implements G4INCL::IAvatar.

Definition at line 1423 of file G4INCLBinaryCollisionAvatar.cc.

1423 {
1424 std::stringstream ss;
1425 ss << "(avatar " << theTime <<" 'nn-collision" << '\n'
1426 << "(list " << '\n'
1427 << particle1->dump()
1428 << particle2->dump()
1429 << "))" << '\n';
1430 return ss.str();
1431 }
G4INCL::Particle::dump
std::string dump() const
Definition G4INCLParticle.hh:1108

◆ getBias()

static G4double G4INCL::BinaryCollisionAvatar::getBias ( )
inlinestatic

Get the global bias factor.

Definition at line 84 of file G4INCLBinaryCollisionAvatar.hh.

84{ return bias; }

◆ getChannel()

G4INCL::IChannel * G4INCL::BinaryCollisionAvatar::getChannel ( )
virtual

Check again the distance of approach. In order for the avatar to be realised, we have to perform a check in the CM system. We define a distance four-vector as

\[ (0, \Delta\vec{x}), \]

where $\Delta\vec{x}$ is the distance vector of the particles at their minimum distance of approach (i.e. at the avatar time). By boosting this four-vector to the CM frame of the two particles and we obtain a new four vector

\[ (\Delta t', \Delta\vec{x}'), \]

with a non-zero time component (the collision happens simultaneously for the two particles in the lab system, but not in the CM system). In order for the avatar to be realised, we require that

\[ |\Delta\vec{x}'| \leq \sqrt{\sigma/\pi}.\]

Note that $|\Delta\vec{x}'|\leq|\Delta\vec{x}|$; thus, the condition above is more restrictive than the check that we perform in G4INCL::Propagation::StandardPropagationModel::generateBinaryCollisionAvatar. In other words, the avatar generation cannot miss any physical collision avatars.

Bias apply for this reaction in order to get the same ParticleBias for all stange particles. Can be reduced after because of the safeguard.

Implements G4INCL::InteractionAvatar.

Definition at line 145 of file G4INCLBinaryCollisionAvatar.cc.

145 {
146 // We already check cutNN at avatar creation time, but we have to check it
147 // again here. For composite projectiles, we might have created independent
148 // avatars with no cutNN before any collision took place.
149 if(particle1->isNucleon()
150 && particle2->isNucleon()
151 && theNucleus->getStore()->getBook().getAcceptedCollisions()!=0) {
152 const G4double energyCM2 = KinematicsUtils::squareTotalEnergyInCM(particle1, particle2);
153 // Below a certain cut value we don't do anything:
154 if(energyCM2 < cutNNSquared) {
155 INCL_DEBUG("CM energy = sqrt(" << energyCM2 << ") MeV < std::sqrt(" << cutNNSquared
156 << ") MeV = cutNN" << "; returning a NULL channel" << '\n');
157 InteractionAvatar::restoreParticles();
158 return NULL;
159 }
160 }
161
162 /** Check again the distance of approach. In order for the avatar to be
163 * realised, we have to perform a check in the CM system. We define a
164 * distance four-vector as
165 * \f[ (0, \Delta\vec{x}), \f]
166 * where \f$\Delta\vec{x}\f$ is the distance vector of the particles at
167 * their minimum distance of approach (i.e. at the avatar time). By
168 * boosting this four-vector to the CM frame of the two particles and we
169 * obtain a new four vector
170 * \f[ (\Delta t', \Delta\vec{x}'), \f]
171 * with a non-zero time component (the collision happens simultaneously for
172 * the two particles in the lab system, but not in the CM system). In order
173 * for the avatar to be realised, we require that
174 * \f[ |\Delta\vec{x}'| \leq \sqrt{\sigma/\pi}.\f]
175 * Note that \f$|\Delta\vec{x}'|\leq|\Delta\vec{x}|\f$; thus, the condition
176 * above is more restrictive than the check that we perform in
177 * G4INCL::Propagation::StandardPropagationModel::generateBinaryCollisionAvatar.
178 * In other words, the avatar generation cannot miss any physical collision
179 * avatars.
180 */
181 ThreeVector minimumDistance = particle1->getPosition();
182 minimumDistance -= particle2->getPosition();
183 const G4double betaDotX = boostVector.dot(minimumDistance);
184 const G4double minDist = Math::tenPi*(minimumDistance.mag2() + betaDotX*betaDotX / (1.-boostVector.mag2()));
185 if(minDist > theCrossSection) {
186 INCL_DEBUG("CM distance of approach is too small: " << minDist << ">" <<
187 theCrossSection <<"; returning a NULL channel" << '\n');
188 InteractionAvatar::restoreParticles();
189 return NULL;
190 }
191
192 /** Bias apply for this reaction in order to get the same
193 * ParticleBias for all stange particles.
194 * Can be reduced after because of the safeguard.
195 */
196 G4double bias_apply = 1.;
197 if(bias != 1.) bias_apply = Particle::getBiasFromVector(Particle::MergeVectorBias(particle1,particle2)) * bias;
198
199//// NN
200 if(particle1->isNucleon() && particle2->isNucleon()) {
201
202 G4double NLKProductionCX = CrossSections::NNToNLK(particle1, particle2)*bias_apply;
203 G4double NSKProductionCX = CrossSections::NNToNSK(particle1, particle2)*bias_apply;
204 G4double NLKpiProductionCX = CrossSections::NNToNLKpi(particle1, particle2)*bias_apply;
205 G4double NSKpiProductionCX = CrossSections::NNToNSKpi(particle1, particle2)*bias_apply;
206 G4double NLK2piProductionCX = CrossSections::NNToNLK2pi(particle1, particle2)*bias_apply;
207 G4double NSK2piProductionCX = CrossSections::NNToNSK2pi(particle1, particle2)*bias_apply;
208 G4double NNKKbProductionCX = CrossSections::NNToNNKKb(particle1, particle2)*bias_apply;
209 G4double NNMissingCX = CrossSections::NNToMissingStrangeness(particle1, particle2)*bias_apply;
210
211 const G4double UnStrangeProdCX = CrossSections::elastic(particle1, particle2) + CrossSections::NNToNDelta(particle1, particle2) + CrossSections::NNToxPiNN(1,particle1, particle2)
212 + CrossSections::NNToxPiNN(2,particle1, particle2) + CrossSections::NNToxPiNN(3,particle1, particle2) + CrossSections::NNToxPiNN(4,particle1, particle2)
213 + CrossSections::NNToNNEtaExclu(particle1, particle2) + CrossSections::NNToNDeltaEta(particle1, particle2) + CrossSections::NNToNNEtaxPi(1,particle1, particle2)
214 + CrossSections::NNToNNEtaxPi(2,particle1, particle2) + CrossSections::NNToNNEtaxPi(3,particle1, particle2) + CrossSections::NNToNNEtaxPi(4,particle1, particle2)
215 + CrossSections::NNToNNOmegaExclu(particle1, particle2) + CrossSections::NNToNDeltaOmega(particle1, particle2) + CrossSections::NNToNNOmegaxPi(1,particle1, particle2)
216 + CrossSections::NNToNNOmegaxPi(2,particle1, particle2) + CrossSections::NNToNNOmegaxPi(3,particle1, particle2) + CrossSections::NNToNNOmegaxPi(4,particle1, particle2);
217 const G4double StrangenessProdCX = (NLKProductionCX + NSKProductionCX + NLKpiProductionCX + NSKpiProductionCX + NLK2piProductionCX + NSK2piProductionCX + NNKKbProductionCX + NNMissingCX)/bias_apply;
218
219 G4double counterweight = (1. - bias_apply * StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX))/(1. - StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX));
220
221 if(counterweight < 0.5) {
222 counterweight = 0.5;
223 bias_apply = 0.5*UnStrangeProdCX/StrangenessProdCX+1;
224 NLKProductionCX = CrossSections::NNToNLK(particle1, particle2)*bias_apply;
225 NSKProductionCX = CrossSections::NNToNSK(particle1, particle2)*bias_apply;
226 NLKpiProductionCX = CrossSections::NNToNLKpi(particle1, particle2)*bias_apply;
227 NSKpiProductionCX = CrossSections::NNToNSKpi(particle1, particle2)*bias_apply;
228 NLK2piProductionCX = CrossSections::NNToNLK2pi(particle1, particle2)*bias_apply;
229 NSK2piProductionCX = CrossSections::NNToNSK2pi(particle1, particle2)*bias_apply;
230 NNKKbProductionCX = CrossSections::NNToNNKKb(particle1, particle2)*bias_apply;
231 NNMissingCX = CrossSections::NNToMissingStrangeness(particle1, particle2)*bias_apply;
232 }
233
234
235 const G4double elasticCX = CrossSections::elastic(particle1, particle2)*counterweight;
236 const G4double deltaProductionCX = CrossSections::NNToNDelta(particle1, particle2)*counterweight;
237 const G4double onePiProductionCX = CrossSections::NNToxPiNN(1,particle1, particle2)*counterweight;
238 const G4double twoPiProductionCX = CrossSections::NNToxPiNN(2,particle1, particle2)*counterweight;
239 const G4double threePiProductionCX = CrossSections::NNToxPiNN(3,particle1, particle2)*counterweight;
240 const G4double fourPiProductionCX = CrossSections::NNToxPiNN(4,particle1, particle2)*counterweight;
241
242 const G4double etaProductionCX = CrossSections::NNToNNEtaExclu(particle1, particle2)*counterweight;
243 const G4double etadeltaProductionCX = CrossSections::NNToNDeltaEta(particle1, particle2)*counterweight;
244 const G4double etaonePiProductionCX = CrossSections::NNToNNEtaxPi(1,particle1, particle2)*counterweight;
245 const G4double etatwoPiProductionCX = CrossSections::NNToNNEtaxPi(2,particle1, particle2)*counterweight;
246 const G4double etathreePiProductionCX = CrossSections::NNToNNEtaxPi(3,particle1, particle2)*counterweight;
247 const G4double etafourPiProductionCX = CrossSections::NNToNNEtaxPi(4,particle1, particle2)*counterweight;
248 const G4double omegaProductionCX = CrossSections::NNToNNOmegaExclu(particle1, particle2)*counterweight;
249 const G4double omegadeltaProductionCX = CrossSections::NNToNDeltaOmega(particle1, particle2)*counterweight;
250 const G4double omegaonePiProductionCX = CrossSections::NNToNNOmegaxPi(1,particle1, particle2)*counterweight;
251 const G4double omegatwoPiProductionCX = CrossSections::NNToNNOmegaxPi(2,particle1, particle2)*counterweight;
252 const G4double omegathreePiProductionCX = CrossSections::NNToNNOmegaxPi(3,particle1, particle2)*counterweight;
253 const G4double omegafourPiProductionCX = CrossSections::NNToNNOmegaxPi(4,particle1, particle2)*counterweight;
254
255 const G4double totCX=CrossSections::total(particle1, particle2);
256
257// assert(std::fabs(totCX-elasticCX-deltaProductionCX-onePiProductionCX-twoPiProductionCX-threePiProductionCX-fourPiProductionCX-NLKProductionCX-NSKProductionCX-NLKpiProductionCX-NSKpiProductionCX-NLK2piProductionCX-NSK2piProductionCX-NNKKbProductionCX-NNMissingCX-etaProductionCX-etadeltaProductionCX-etaonePiProductionCX-etatwoPiProductionCX-etathreePiProductionCX-etafourPiProductionCX-omegaProductionCX-omegadeltaProductionCX-omegaonePiProductionCX-omegatwoPiProductionCX-omegathreePiProductionCX-omegafourPiProductionCX) < 0.5);
258
259 const G4double rChannel=Random::shoot() * totCX;
260
261 if(elasticCX > rChannel) {
262// Elastic NN channel
263 isElastic = true;
264 INCL_DEBUG("NN interaction: elastic channel chosen" << '\n');
265 weight = counterweight;
266 return new ElasticChannel(particle1, particle2);
267 } else if((elasticCX + deltaProductionCX) > rChannel) {
268 isElastic = false;
269// NN -> N Delta channel is chosen
270 INCL_DEBUG("NN interaction: Delta channel chosen" << '\n');
271 weight = counterweight;
272 return new DeltaProductionChannel(particle1, particle2);
273 } else if(elasticCX + deltaProductionCX + onePiProductionCX > rChannel) {
274 isElastic = false;
275// NN -> PiNN channel is chosen
276 INCL_DEBUG("NN interaction: one Pion channel chosen" << '\n');
277 weight = counterweight;
278 return new NNToMultiPionsChannel(1,particle1, particle2);
279 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX > rChannel) {
280 isElastic = false;
281// NN -> 2PiNN channel is chosen
282 INCL_DEBUG("NN interaction: two Pions channel chosen" << '\n');
283 weight = counterweight;
284 return new NNToMultiPionsChannel(2,particle1, particle2);
285 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX > rChannel) {
286 isElastic = false;
287// NN -> 3PiNN channel is chosen
288 INCL_DEBUG("NN interaction: three Pions channel chosen" << '\n');
289 weight = counterweight;
290 return new NNToMultiPionsChannel(3,particle1, particle2);
291 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX > rChannel) {
292 isElastic = false;
293// NN -> 4PiNN channel is chosen
294 INCL_DEBUG("NN interaction: four Pions channel chosen" << '\n');
295 weight = counterweight;
296 return new NNToMultiPionsChannel(4,particle1, particle2);
297 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
298 + etaProductionCX > rChannel) {
299 isElastic = false;
300// NN -> NNEta channel is chosen
301 INCL_DEBUG("NN interaction: Eta channel chosen" << '\n');
302 weight = counterweight;
303 return new NNToNNEtaChannel(particle1, particle2);
304 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
305 + etaProductionCX + etadeltaProductionCX > rChannel) {
306 isElastic = false;
307// NN -> N Delta Eta channel is chosen
308 INCL_DEBUG("NN interaction: Delta Eta channel chosen" << '\n');
309 weight = counterweight;
310 return new NDeltaEtaProductionChannel(particle1, particle2);
311 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
312 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX > rChannel) {
313 isElastic = false;
314// NN -> EtaPiNN channel is chosen
315 INCL_DEBUG("NN interaction: Eta + one Pion channel chosen" << '\n');
316 weight = counterweight;
317 return new NNEtaToMultiPionsChannel(1,particle1, particle2);
318 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
319 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX > rChannel) {
320 isElastic = false;
321// NN -> Eta2PiNN channel is chosen
322 INCL_DEBUG("NN interaction: Eta + two Pions channel chosen" << '\n');
323 weight = counterweight;
324 return new NNEtaToMultiPionsChannel(2,particle1, particle2);
325 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
326 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX > rChannel) {
327 isElastic = false;
328// NN -> Eta3PiNN channel is chosen
329 INCL_DEBUG("NN interaction: Eta + three Pions channel chosen" << '\n');
330 weight = counterweight;
331 return new NNEtaToMultiPionsChannel(3,particle1, particle2);
332 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
333 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX > rChannel) {
334 isElastic = false;
335// NN -> Eta4PiNN channel is chosen
336 INCL_DEBUG("NN interaction: Eta + four Pions channel chosen" << '\n');
337 weight = counterweight;
338 return new NNEtaToMultiPionsChannel(4,particle1, particle2);
339 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
340 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
341 + omegaProductionCX > rChannel) {
342 isElastic = false;
343// NN -> NNOmega channel is chosen
344 INCL_DEBUG("NN interaction: Omega channel chosen" << '\n');
345 weight = counterweight;
346 return new NNToNNOmegaChannel(particle1, particle2);
347 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
348 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
349 + omegaProductionCX + omegadeltaProductionCX > rChannel) {
350 isElastic = false;
351// NN -> N Delta Omega channel is chosen
352 INCL_DEBUG("NN interaction: Delta Omega channel chosen" << '\n');
353 weight = counterweight;
354 return new NDeltaOmegaProductionChannel(particle1, particle2);
355 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
356 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
357 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX > rChannel) {
358 isElastic = false;
359// NN -> OmegaPiNN channel is chosen
360 INCL_DEBUG("NN interaction: Omega + one Pion channel chosen" << '\n');
361 weight = counterweight;
362 return new NNOmegaToMultiPionsChannel(1,particle1, particle2);
363 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
364 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
365 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX > rChannel) {
366 isElastic = false;
367// NN -> Omega2PiNN channel is chosen
368 INCL_DEBUG("NN interaction: Omega + two Pions channel chosen" << '\n');
369 weight = counterweight;
370 return new NNOmegaToMultiPionsChannel(2,particle1, particle2);
371 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
372 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
373 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX > rChannel) {
374 isElastic = false;
375// NN -> Omega3PiNN channel is chosen
376 INCL_DEBUG("NN interaction: Omega + three Pions channel chosen" << '\n');
377 weight = counterweight;
378 return new NNOmegaToMultiPionsChannel(3,particle1, particle2);
379 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
380 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
381 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX > rChannel) {
382 isElastic = false;
383// NN -> Omega4PiNN channel is chosen
384 INCL_DEBUG("NN interaction: Omega + four Pions channel chosen" << '\n');
385 weight = counterweight;
386 return new NNOmegaToMultiPionsChannel(4,particle1, particle2);
387 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
388 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
389 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
390 + NLKProductionCX > rChannel) {
391 isElastic = false;
392 isStrangeProduction = true;
393// NN -> NLK channel is chosen
394 INCL_DEBUG("NN interaction: NLK channel chosen" << '\n');
395 weight = bias_apply;
396 return new NNToNLKChannel(particle1, particle2);
397 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
398 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
399 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
400 + NLKProductionCX + NLKpiProductionCX > rChannel) {
401 isElastic = false;
402 isStrangeProduction = true;
403// NN -> NLKpi channel is chosen
404 INCL_DEBUG("NN interaction: NLKpi channel chosen" << '\n');
405 weight = bias_apply;
406 return new NNToNLKpiChannel(particle1, particle2);
407 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
408 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
409 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
410 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX > rChannel) {
411 isElastic = false;
412 isStrangeProduction = true;
413// NN -> NLK2pi channel is chosen
414 INCL_DEBUG("NN interaction: NLK2pi channel chosen" << '\n');
415 weight = bias_apply;
416 return new NNToNLK2piChannel(particle1, particle2);
417 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
418 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
419 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
420 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX + NSKProductionCX > rChannel) {
421 isElastic = false;
422 isStrangeProduction = true;
423// NN -> NSK channel is chosen
424 INCL_DEBUG("NN interaction: NSK channel chosen" << '\n');
425 weight = bias_apply;
426 return new NNToNSKChannel(particle1, particle2);
427 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
428 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
429 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
430 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX + NSKProductionCX + NSKpiProductionCX > rChannel) {
431 isElastic = false;
432 isStrangeProduction = true;
433// NN -> NSKpi channel is chosen
434 INCL_DEBUG("NN interaction: NSKpi channel chosen" << '\n');
435 weight = bias_apply;
436 return new NNToNSKpiChannel(particle1, particle2);
437 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
438 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
439 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
440 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX + NSKProductionCX + NSKpiProductionCX + NSK2piProductionCX > rChannel) {
441 isElastic = false;
442 isStrangeProduction = true;
443// NN -> NSK2pi channel is chosen
444 INCL_DEBUG("NN interaction: NSK2pi channel chosen" << '\n');
445 weight = bias_apply;
446 return new NNToNSK2piChannel(particle1, particle2);
447 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
448 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
449 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
450 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX + NSKProductionCX + NSKpiProductionCX + NSK2piProductionCX + NNKKbProductionCX > rChannel) {
451 isElastic = false;
452 isStrangeProduction = true;
453// NN -> NNKKb channel is chosen
454 INCL_DEBUG("NN interaction: NNKKb channel chosen" << '\n');
455 weight = bias_apply;
456 return new NNToNNKKbChannel(particle1, particle2);
457 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + fourPiProductionCX
458 + etaProductionCX + etadeltaProductionCX + etaonePiProductionCX + etatwoPiProductionCX + etathreePiProductionCX + etafourPiProductionCX
459 + omegaProductionCX + omegadeltaProductionCX + omegaonePiProductionCX + omegatwoPiProductionCX + omegathreePiProductionCX + omegafourPiProductionCX
460 + NLKProductionCX + NLKpiProductionCX + NLK2piProductionCX + NSKProductionCX + NSKpiProductionCX + NSK2piProductionCX + NNKKbProductionCX + NNMissingCX> rChannel) {
461 isElastic = false;
462 isStrangeProduction = true;
463// NN -> Missing Strangeness channel is chosen
464 INCL_DEBUG("NN interaction: Missing Strangeness channel chosen" << '\n');
465 weight = bias_apply;
466 return new NNToMissingStrangenessChannel(particle1, particle2);
467 } else {
468 INCL_WARN("inconsistency within the NN Cross Sections (sum!=inelastic)" << '\n');
469 if(NNMissingCX>0.) {
470 INCL_WARN("Returning an Missing Strangeness channel" << '\n');
471 weight = bias_apply;
472 isElastic = false;
473 isStrangeProduction = true;
474 return new NNToNNKKbChannel(particle1, particle2);
475 } else if(NNKKbProductionCX>0.) {
476 INCL_WARN("Returning an NNKKb channel" << '\n');
477 weight = bias_apply;
478 isElastic = false;
479 isStrangeProduction = true;
480 return new NNToNNKKbChannel(particle1, particle2);
481 } else if(NSK2piProductionCX>0.) {
482 INCL_WARN("Returning an NSK2pi channel" << '\n');
483 weight = bias_apply;
484 isElastic = false;
485 isStrangeProduction = true;
486 return new NNToNSK2piChannel(particle1, particle2);
487 } else if(NSKpiProductionCX>0.) {
488 INCL_WARN("Returning an NSKpi channel" << '\n');
489 weight = bias_apply;
490 isElastic = false;
491 isStrangeProduction = true;
492 return new NNToNSKpiChannel(particle1, particle2);
493 } else if(NSKProductionCX>0.) {
494 INCL_WARN("Returning an NSK channel" << '\n');
495 weight = bias_apply;
496 isElastic = false;
497 isStrangeProduction = true;
498 return new NNToNSKChannel(particle1, particle2);
499 } else if(NLK2piProductionCX>0.) {
500 INCL_WARN("Returning an NLK2pi channel" << '\n');
501 weight = bias_apply;
502 isElastic = false;
503 isStrangeProduction = true;
504 return new NNToNLK2piChannel(particle1, particle2);
505 } else if(NLKpiProductionCX>0.) {
506 INCL_WARN("Returning an NLKpi channel" << '\n');
507 weight = bias_apply;
508 isElastic = false;
509 isStrangeProduction = true;
510 return new NNToNLKpiChannel(particle1, particle2);
511 } else if(NLKProductionCX>0.) {
512 INCL_WARN("Returning an NLK channel" << '\n');
513 weight = bias_apply;
514 isElastic = false;
515 isStrangeProduction = true;
516 return new NNToNLKChannel(particle1, particle2);
517 } else if(omegafourPiProductionCX>0.) {
518 INCL_WARN("Returning an Omega + four Pions channel" << '\n');
519 weight = counterweight;
520 isElastic = false;
521 return new NNOmegaToMultiPionsChannel(4,particle1, particle2);
522 } else if(omegathreePiProductionCX>0.) {
523 INCL_WARN("Returning an Omega + three Pions channel" << '\n');
524 weight = counterweight;
525 isElastic = false;
526 return new NNOmegaToMultiPionsChannel(3,particle1, particle2);
527 } else if(omegatwoPiProductionCX>0.) {
528 INCL_WARN("Returning an Omega + two Pions channel" << '\n');
529 weight = counterweight;
530 isElastic = false;
531 return new NNOmegaToMultiPionsChannel(2,particle1, particle2);
532 } else if(omegaonePiProductionCX>0.) {
533 INCL_WARN("Returning an Omega + one Pion channel" << '\n');
534 weight = counterweight;
535 isElastic = false;
536 return new NNOmegaToMultiPionsChannel(1,particle1, particle2);
537 } else if(omegadeltaProductionCX>0.) {
538 INCL_WARN("Returning an Omega + Delta channel" << '\n');
539 weight = counterweight;
540 isElastic = false;
541 return new NDeltaOmegaProductionChannel(particle1, particle2);
542 } else if(omegaProductionCX>0.) {
543 INCL_WARN("Returning an Omega channel" << '\n');
544 weight = counterweight;
545 isElastic = false;
546 return new NNToNNOmegaChannel(particle1, particle2);
547 } else if(etafourPiProductionCX>0.) {
548 INCL_WARN("Returning an Eta + four Pions channel" << '\n');
549 weight = counterweight;
550 isElastic = false;
551 return new NNEtaToMultiPionsChannel(4,particle1, particle2);
552 } else if(etathreePiProductionCX>0.) {
553 INCL_WARN("Returning an Eta + threev channel" << '\n');
554 weight = counterweight;
555 isElastic = false;
556 return new NNEtaToMultiPionsChannel(3,particle1, particle2);
557 } else if(etatwoPiProductionCX>0.) {
558 INCL_WARN("Returning an Eta + two Pions channel" << '\n');
559 weight = counterweight;
560 isElastic = false;
561 return new NNEtaToMultiPionsChannel(2,particle1, particle2);
562 } else if(etaonePiProductionCX>0.) {
563 INCL_WARN("Returning an Eta + one Pion channel" << '\n');
564 weight = counterweight;
565 isElastic = false;
566 return new NNEtaToMultiPionsChannel(1,particle1, particle2);
567 } else if(etadeltaProductionCX>0.) {
568 INCL_WARN("Returning an Eta + Delta channel" << '\n');
569 weight = counterweight;
570 isElastic = false;
571 return new NDeltaEtaProductionChannel(particle1, particle2);
572 } else if(etaProductionCX>0.) {
573 INCL_WARN("Returning an Eta channel" << '\n');
574 weight = counterweight;
575 isElastic = false;
576 return new NNToNNEtaChannel(particle1, particle2);
577 } else if(fourPiProductionCX>0.) {
578 INCL_WARN("Returning a 4pi channel" << '\n');
579 weight = counterweight;
580 isElastic = false;
581 return new NNToMultiPionsChannel(4,particle1, particle2);
582 } else if(threePiProductionCX>0.) {
583 INCL_WARN("Returning a 3pi channel" << '\n');
584 weight = counterweight;
585 isElastic = false;
586 return new NNToMultiPionsChannel(3,particle1, particle2);
587 } else if(twoPiProductionCX>0.) {
588 INCL_WARN("Returning a 2pi channel" << '\n');
589 weight = counterweight;
590 isElastic = false;
591 return new NNToMultiPionsChannel(2,particle1, particle2);
592 } else if(onePiProductionCX>0.) {
593 INCL_WARN("Returning a 1pi channel" << '\n');
594 weight = counterweight;
595 isElastic = false;
596 return new NNToMultiPionsChannel(1,particle1, particle2);
597 } else if(deltaProductionCX>0.) {
598 INCL_WARN("Returning a delta-production channel" << '\n');
599 weight = counterweight;
600 isElastic = false;
601 return new DeltaProductionChannel(particle1, particle2);
602 } else {
603 INCL_WARN("Returning an elastic channel" << '\n');
604 weight = counterweight;
605 isElastic = true;
606 return new ElasticChannel(particle1, particle2);
607 }
608 }
609
610//// NDelta
611 }
612 else if((particle1->isNucleon() && particle2->isDelta()) ||
613 (particle1->isDelta() && particle2->isNucleon())) {
614
615 G4double NLKProductionCX = CrossSections::NDeltaToNLK(particle1, particle2)*bias_apply;
616 G4double NSKProductionCX = CrossSections::NDeltaToNSK(particle1, particle2)*bias_apply;
617 G4double DeltaLKProductionCX = CrossSections::NDeltaToDeltaLK(particle1, particle2)*bias_apply;
618 G4double DeltaSKProductionCX = CrossSections::NDeltaToDeltaSK(particle1, particle2)*bias_apply;
619 G4double NNKKbProductionCX = CrossSections::NDeltaToNNKKb(particle1, particle2)*bias_apply;
620
621 const G4double UnStrangeProdCX = CrossSections::elastic(particle1, particle2) + CrossSections::NDeltaToNN(particle1, particle2);
622 const G4double StrangenessProdCX = (NLKProductionCX + NSKProductionCX + DeltaLKProductionCX + DeltaSKProductionCX + NNKKbProductionCX)/bias_apply;
623
624 G4double counterweight = (1. - bias_apply * StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX))/(1. - StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX));
625
626 if(counterweight < 0.5){
627 counterweight = 0.5;
628 bias_apply = 0.5*UnStrangeProdCX/StrangenessProdCX+1;
629
630 NLKProductionCX = CrossSections::NDeltaToNLK(particle1, particle2)*bias_apply;
631 NSKProductionCX = CrossSections::NDeltaToNSK(particle1, particle2)*bias_apply;
632 DeltaLKProductionCX = CrossSections::NDeltaToDeltaLK(particle1, particle2)*bias_apply;
633 DeltaSKProductionCX = CrossSections::NDeltaToDeltaSK(particle1, particle2)*bias_apply;
634 NNKKbProductionCX = CrossSections::NDeltaToNNKKb(particle1, particle2)*bias_apply;
635 }
636
637 G4double elasticCX = CrossSections::elastic(particle1, particle2)*counterweight;
638 G4double recombinationCX = CrossSections::NDeltaToNN(particle1, particle2)*counterweight;
639
640 const G4double rChannel=Random::shoot() * (StrangenessProdCX + UnStrangeProdCX);
641
642 if(elasticCX > rChannel) {
643 isElastic = true;
644// Elastic N Delta channel
645 INCL_DEBUG("NDelta interaction: elastic channel chosen" << '\n');
646 weight = counterweight;
647 return new ElasticChannel(particle1, particle2);
648 } else if (elasticCX + recombinationCX > rChannel){
649 isElastic = false;
650// Recombination
651// NDelta -> NN channel is chosen
652 INCL_DEBUG("NDelta interaction: recombination channel chosen" << '\n');
653 weight = counterweight;
654 return new RecombinationChannel(particle1, particle2);
655 } else if (elasticCX + recombinationCX + NLKProductionCX > rChannel){
656 isElastic = false;
657 isStrangeProduction = true;
658// NDelta -> NLK channel is chosen
659 INCL_DEBUG("NDelta interaction: NLK channel chosen" << '\n');
660 weight = bias_apply;
661 return new NDeltaToNLKChannel(particle1, particle2);
662 } else if (elasticCX + recombinationCX + NLKProductionCX + NSKProductionCX > rChannel){
663 isElastic = false;
664 isStrangeProduction = true;
665// NDelta -> NSK channel is chosen
666 INCL_DEBUG("NDelta interaction: NSK channel chosen" << '\n');
667 weight = bias_apply;
668 return new NDeltaToNSKChannel(particle1, particle2);
669 } else if (elasticCX + recombinationCX + NLKProductionCX + NSKProductionCX + DeltaLKProductionCX > rChannel){
670 isElastic = false;
671 isStrangeProduction = true;
672// NDelta -> DeltaLK channel is chosen
673 INCL_DEBUG("NDelta interaction: DeltaLK channel chosen" << '\n');
674 weight = bias_apply;
675 return new NDeltaToDeltaLKChannel(particle1, particle2);
676 } else if (elasticCX + recombinationCX + NLKProductionCX + NSKProductionCX + DeltaLKProductionCX + DeltaSKProductionCX > rChannel){
677 isElastic = false;
678 isStrangeProduction = true;
679// NDelta -> DeltaSK channel is chosen
680 INCL_DEBUG("NDelta interaction: DeltaSK channel chosen" << '\n');
681 weight = bias_apply;
682 return new NDeltaToDeltaSKChannel(particle1, particle2);
683 } else if (elasticCX + recombinationCX + NLKProductionCX + NSKProductionCX + DeltaLKProductionCX + DeltaSKProductionCX + NNKKbProductionCX > rChannel){
684 isElastic = false;
685 isStrangeProduction = true;
686// NDelta -> NNKKb channel is chosen
687 INCL_DEBUG("NDelta interaction: NNKKb channel chosen" << '\n');
688 weight = bias_apply;
689 return new NDeltaToNNKKbChannel(particle1, particle2);
690 }
691 else{
692 INCL_ERROR("rChannel > (StrangenessProdCX + UnStrangeProdCX) in NDelta interaction: return an elastic channel" << '\n');
693 weight = counterweight;
694 isElastic = true;
695 return new ElasticChannel(particle1, particle2);
696 }
697
698//// DeltaDelta
699 } else if(particle1->isDelta() && particle2->isDelta()) {
700 isElastic = true;
701 INCL_DEBUG("DeltaDelta interaction: elastic channel chosen" << '\n');
702 return new ElasticChannel(particle1, particle2);
703
704//// PiN
705 } else if(isPiN) {
706
707 G4double LKProdCX = CrossSections::NpiToLK(particle1,particle2)*bias_apply;
708 G4double SKProdCX = CrossSections::NpiToSK(particle1,particle2)*bias_apply;
709 G4double LKpiProdCX = CrossSections::NpiToLKpi(particle1,particle2)*bias_apply;
710 G4double SKpiProdCX = CrossSections::NpiToSKpi(particle1,particle2)*bias_apply;
711 G4double LK2piProdCX = CrossSections::NpiToLK2pi(particle1,particle2)*bias_apply;
712 G4double SK2piProdCX = CrossSections::NpiToSK2pi(particle1,particle2)*bias_apply;
713 G4double NKKbProdCX = CrossSections::NpiToNKKb(particle1,particle2)*bias_apply;
714 G4double MissingCX = CrossSections::NpiToMissingStrangeness(particle1,particle2)*bias_apply;
715
716 const G4double UnStrangeProdCX = CrossSections::elastic(particle1, particle2) + CrossSections::piNToDelta(particle1, particle2)
717 + CrossSections::piNToxPiN(2,particle1, particle2) + CrossSections::piNToxPiN(3,particle1, particle2) + CrossSections::piNToxPiN(4,particle1, particle2)
718 + CrossSections::piNToEtaN(particle1, particle2) + CrossSections::piNToOmegaN(particle1, particle2);
719 const G4double StrangenessProdCX = (LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX + LK2piProdCX + SK2piProdCX + NKKbProdCX + MissingCX)/bias_apply;
720
721 G4double counterweight = (1. - bias_apply * StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX))/(1. - StrangenessProdCX / (StrangenessProdCX + UnStrangeProdCX));
722
723 if(counterweight < 0.5) {
724 counterweight = 0.5;
725 bias_apply = 0.5*UnStrangeProdCX/StrangenessProdCX+1;
726 LKProdCX = CrossSections::NpiToLK(particle1,particle2)*bias_apply;
727 SKProdCX = CrossSections::NpiToSK(particle1,particle2)*bias_apply;
728 LKpiProdCX = CrossSections::NpiToLKpi(particle1,particle2)*bias_apply;
729 SKpiProdCX = CrossSections::NpiToSKpi(particle1,particle2)*bias_apply;
730 LK2piProdCX = CrossSections::NpiToLK2pi(particle1,particle2)*bias_apply;
731 SK2piProdCX = CrossSections::NpiToSK2pi(particle1,particle2)*bias_apply;
732 NKKbProdCX = CrossSections::NpiToNKKb(particle1,particle2)*bias_apply;
733 MissingCX = CrossSections::NpiToMissingStrangeness(particle1,particle2)*bias_apply;
734 }
735
736
737 const G4double elasticCX = CrossSections::elastic(particle1, particle2)*counterweight;
738 const G4double deltaProductionCX = CrossSections::piNToDelta(particle1, particle2)*counterweight;
739 const G4double onePiProductionCX = CrossSections::piNToxPiN(2,particle1, particle2)*counterweight;
740 const G4double twoPiProductionCX = CrossSections::piNToxPiN(3,particle1, particle2)*counterweight;
741 const G4double threePiProductionCX = CrossSections::piNToxPiN(4,particle1, particle2)*counterweight;
742 const G4double etaProductionCX = CrossSections::piNToEtaN(particle1, particle2)*counterweight;
743 const G4double omegaProductionCX = CrossSections::piNToOmegaN(particle1, particle2)*counterweight;
744
745 const G4double totCX=CrossSections::total(particle1, particle2);
746
747// assert(std::fabs(totCX-elasticCX-deltaProductionCX-onePiProductionCX-twoPiProductionCX-threePiProductionCX-etaProductionCX-omegaProductionCX-LKProdCX-SKProdCX-LKpiProdCX-SKpiProdCX-LK2piProdCX-SK2piProdCX-NKKbProdCX-MissingCX) < 0.15);
748
749 const G4double rChannel=Random::shoot() * totCX;
750
751 if(elasticCX > rChannel) {
752 isElastic = true;
753// Elastic PiN channel
754 INCL_DEBUG("PiN interaction: elastic channel chosen" << '\n');
755 weight = counterweight;
756 return new PiNElasticChannel(particle1, particle2);
757 } else if(elasticCX + deltaProductionCX > rChannel) {
758 isElastic = false;
759// PiN -> Delta channel is chosen
760 INCL_DEBUG("PiN interaction: Delta channel chosen" << '\n');
761 weight = counterweight;
762 return new PiNToDeltaChannel(particle1, particle2);
763 } else if(elasticCX + deltaProductionCX + onePiProductionCX > rChannel) {
764 isElastic = false;
765// PiN -> PiNPi channel is chosen
766 INCL_DEBUG("PiN interaction: one Pion channel chosen" << '\n');
767 weight = counterweight;
768 return new PiNToMultiPionsChannel(2,particle1, particle2);
769 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX > rChannel) {
770 isElastic = false;
771// PiN -> PiN2Pi channel is chosen
772 INCL_DEBUG("PiN interaction: two Pions channel chosen" << '\n');
773 weight = counterweight;
774 return new PiNToMultiPionsChannel(3,particle1, particle2);
775 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX > rChannel) {
776 isElastic = false;
777// PiN -> PiN3Pi channel is chosen
778 INCL_DEBUG("PiN interaction: three Pions channel chosen" << '\n');
779 weight = counterweight;
780 return new PiNToMultiPionsChannel(4,particle1, particle2);
781 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX > rChannel) {
782 isElastic = false;
783// PiN -> EtaN channel is chosen
784 INCL_DEBUG("PiN interaction: Eta channel chosen" << '\n');
785 weight = counterweight;
786 return new PiNToEtaChannel(particle1, particle2);
787 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX > rChannel) {
788 isElastic = false;
789// PiN -> OmegaN channel is chosen
790 INCL_DEBUG("PiN interaction: Omega channel chosen" << '\n');
791 weight = counterweight;
792 return new PiNToOmegaChannel(particle1, particle2);
793 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
794 + LKProdCX > rChannel) {
795 isElastic = false;
796 isStrangeProduction = true;
797// PiN -> LK channel is chosen
798 INCL_DEBUG("PiN interaction: LK channel chosen" << '\n');
799 weight = bias_apply;
800 return new NpiToLKChannel(particle1, particle2);
801 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
802 + LKProdCX + SKProdCX > rChannel) {
803 isElastic = false;
804 isStrangeProduction = true;
805// PiN -> SK channel is chosen
806 INCL_DEBUG("PiN interaction: SK channel chosen" << '\n');
807 weight = bias_apply;
808 return new NpiToSKChannel(particle1, particle2);
809 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
810 + LKProdCX + SKProdCX + LKpiProdCX > rChannel) {
811 isElastic = false;
812 isStrangeProduction = true;
813// PiN -> LKpi channel is chosen
814 INCL_DEBUG("PiN interaction: LKpi channel chosen" << '\n');
815 weight = bias_apply;
816 return new NpiToLKpiChannel(particle1, particle2);
817 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
818 + LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX > rChannel) {
819 isElastic = false;
820 isStrangeProduction = true;
821// PiN -> SKpi channel is chosen
822 INCL_DEBUG("PiN interaction: SKpi channel chosen" << '\n');
823 weight = bias_apply;
824 return new NpiToSKpiChannel(particle1, particle2);
825 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
826 + LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX + LK2piProdCX > rChannel) {
827 isElastic = false;
828 isStrangeProduction = true;
829// PiN -> LK2pi channel is chosen
830 INCL_DEBUG("PiN interaction: LK2pi channel chosen" << '\n');
831 weight = bias_apply;
832 return new NpiToLK2piChannel(particle1, particle2);
833 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
834 + LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX + LK2piProdCX + SK2piProdCX > rChannel) {
835 isElastic = false;
836 isStrangeProduction = true;
837// PiN -> SK2pi channel is chosen
838 INCL_DEBUG("PiN interaction: SK2pi channel chosen" << '\n');
839 weight = bias_apply;
840 return new NpiToSK2piChannel(particle1, particle2);
841 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
842 + LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX + LK2piProdCX + SK2piProdCX + NKKbProdCX > rChannel) {
843 isElastic = false;
844 isStrangeProduction = true;
845// PiN -> NKKb channel is chosen
846 INCL_DEBUG("PiN interaction: NKKb channel chosen" << '\n');
847 weight = bias_apply;
848 return new NpiToNKKbChannel(particle1, particle2);
849 } else if(elasticCX + deltaProductionCX + onePiProductionCX + twoPiProductionCX + threePiProductionCX + etaProductionCX+ omegaProductionCX
850 + LKProdCX + SKProdCX + LKpiProdCX + SKpiProdCX + LK2piProdCX + SK2piProdCX + NKKbProdCX + MissingCX> rChannel) {
851 isElastic = false;
852 isStrangeProduction = true;
853// PiN -> Missinge Strangeness channel is chosen
854 INCL_DEBUG("PiN interaction: Missinge Strangeness channel chosen" << '\n');
855 weight = bias_apply;
856 return new NpiToMissingStrangenessChannel(particle1, particle2);
857 }
858 else {
859 INCL_WARN("inconsistency within the PiN Cross Sections (sum!=inelastic)" << '\n');
860 if(MissingCX>0.) {
861 INCL_WARN("Returning a Missinge Strangeness channel" << '\n');
862 weight = bias_apply;
863 isElastic = false;
864 isStrangeProduction = true;
865 return new NpiToMissingStrangenessChannel(particle1, particle2);
866 } else if(NKKbProdCX>0.) {
867 INCL_WARN("Returning a NKKb channel" << '\n');
868 weight = bias_apply;
869 isElastic = false;
870 isStrangeProduction = true;
871 return new NpiToNKKbChannel(particle1, particle2);
872 } else if(SK2piProdCX>0.) {
873 INCL_WARN("Returning a SK2pi channel" << '\n');
874 weight = bias_apply;
875 isElastic = false;
876 isStrangeProduction = true;
877 return new NpiToSK2piChannel(particle1, particle2);
878 } else if(LK2piProdCX>0.) {
879 INCL_WARN("Returning a LK2pi channel" << '\n');
880 weight = bias_apply;
881 isElastic = false;
882 isStrangeProduction = true;
883 return new NpiToLK2piChannel(particle1, particle2);
884 } else if(SKpiProdCX>0.) {
885 INCL_WARN("Returning a SKpi channel" << '\n');
886 weight = bias_apply;
887 isElastic = false;
888 isStrangeProduction = true;
889 return new NpiToSKpiChannel(particle1, particle2);
890 } else if(LKpiProdCX>0.) {
891 INCL_WARN("Returning a LKpi channel" << '\n');
892 weight = bias_apply;
893 isElastic = false;
894 isStrangeProduction = true;
895 return new NpiToLKpiChannel(particle1, particle2);
896 } else if(SKProdCX>0.) {
897 INCL_WARN("Returning a SK channel" << '\n');
898 weight = bias_apply;
899 isElastic = false;
900 isStrangeProduction = true;
901 return new NpiToSKChannel(particle1, particle2);
902 } else if(LKProdCX>0.) {
903 INCL_WARN("Returning a LK channel" << '\n');
904 weight = bias_apply;
905 isElastic = false;
906 isStrangeProduction = true;
907 return new NpiToLKChannel(particle1, particle2);
908 } else if(omegaProductionCX>0.) {
909 INCL_WARN("Returning a Omega channel" << '\n');
910 weight = counterweight;
911 isElastic = false;
912 return new PiNToOmegaChannel(particle1, particle2);
913 } else if(etaProductionCX>0.) {
914 INCL_WARN("Returning a Eta channel" << '\n');
915 weight = counterweight;
916 isElastic = false;
917 return new PiNToEtaChannel(particle1, particle2);
918 } else if(threePiProductionCX>0.) {
919 INCL_WARN("Returning a 3pi channel" << '\n');
920 weight = counterweight;
921 isElastic = false;
922 return new PiNToMultiPionsChannel(4,particle1, particle2);
923 } else if(twoPiProductionCX>0.) {
924 INCL_WARN("Returning a 2pi channel" << '\n');
925 weight = counterweight;
926 isElastic = false;
927 return new PiNToMultiPionsChannel(3,particle1, particle2);
928 } else if(onePiProductionCX>0.) {
929 INCL_WARN("Returning a 1pi channel" << '\n');
930 weight = counterweight;
931 isElastic = false;
932 return new PiNToMultiPionsChannel(2,particle1, particle2);
933 } else if(deltaProductionCX>0.) {
934 INCL_WARN("Returning a delta-production channel" << '\n');
935 weight = counterweight;
936 isElastic = false;
937 return new PiNToDeltaChannel(particle1, particle2);
938 } else {
939 INCL_WARN("Returning an elastic channel" << '\n');
940 weight = counterweight;
941 isElastic = true;
942 return new PiNElasticChannel(particle1, particle2);
943 }
944 }
945 } else if ((particle1->isNucleon() && particle2->isEta()) || (particle2->isNucleon() && particle1->isEta())) {
946//// EtaN
947
948 const G4double elasticCX = CrossSections::elastic(particle1, particle2);
949 const G4double onePiProductionCX = CrossSections::etaNToPiN(particle1, particle2);
950 const G4double twoPiProductionCX = CrossSections::etaNToPiPiN(particle1, particle2);
951 const G4double totCX=CrossSections::total(particle1, particle2);
952// assert(std::fabs(totCX-elasticCX-onePiProductionCX-twoPiProductionCX)<1.);
953
954 const G4double rChannel=Random::shoot() * totCX;
955
956 if(elasticCX > rChannel) {
957// Elastic EtaN channel
958 isElastic = true;
959 INCL_DEBUG("EtaN interaction: elastic channel chosen" << '\n');
960 return new EtaNElasticChannel(particle1, particle2);
961 } else if(elasticCX + onePiProductionCX > rChannel) {
962 isElastic = false;
963// EtaN -> EtaPiN channel is chosen
964 INCL_DEBUG("EtaN interaction: PiN channel chosen" << '\n');
965 return new EtaNToPiNChannel(particle1, particle2);
966 } else if(elasticCX + onePiProductionCX + twoPiProductionCX > rChannel) {
967 isElastic = false;
968// EtaN -> EtaPiPiN channel is chosen
969 INCL_DEBUG("EtaN interaction: PiPiN channel chosen" << '\n');
970 return new EtaNToPiPiNChannel(particle1, particle2);
971 }
972
973 else {
974 INCL_WARN("inconsistency within the EtaN Cross Sections (sum!=inelastic)" << '\n');
975 if(twoPiProductionCX>0.) {
976 INCL_WARN("Returning a PiPiN channel" << '\n');
977 isElastic = false;
978 return new EtaNToPiPiNChannel(particle1, particle2);
979 } else if(onePiProductionCX>0.) {
980 INCL_WARN("Returning a PiN channel" << '\n');
981 isElastic = false;
982 return new EtaNToPiNChannel(particle1, particle2);
983 } else {
984 INCL_WARN("Returning an elastic channel" << '\n');
985 isElastic = true;
986 return new EtaNElasticChannel(particle1, particle2);
987 }
988 }
989
990 } else if ((particle1->isNucleon() && particle2->isOmega()) || (particle2->isNucleon() && particle1->isOmega())) {
991//// OmegaN
992
993 const G4double elasticCX = CrossSections::elastic(particle1, particle2);
994 const G4double onePiProductionCX = CrossSections::omegaNToPiN(particle1, particle2);
995 const G4double twoPiProductionCX = CrossSections::omegaNToPiPiN(particle1, particle2);
996 const G4double totCX=CrossSections::total(particle1, particle2);
997// assert(std::fabs(totCX-elasticCX-onePiProductionCX-twoPiProductionCX)<1.);
998
999 const G4double rChannel=Random::shoot() * totCX;
1000
1001 if(elasticCX > rChannel) {
1002// Elastic OmegaN channel
1003 isElastic = true;
1004 INCL_DEBUG("OmegaN interaction: elastic channel chosen" << '\n');
1005 return new OmegaNElasticChannel(particle1, particle2);
1006 } else if(elasticCX + onePiProductionCX > rChannel) {
1007 isElastic = false;
1008// OmegaN -> PiN channel is chosen
1009 INCL_DEBUG("OmegaN interaction: PiN channel chosen" << '\n');
1010 return new OmegaNToPiNChannel(particle1, particle2);
1011 } else if(elasticCX + onePiProductionCX + twoPiProductionCX > rChannel) {
1012 isElastic = false;
1013// OmegaN -> PiPiN channel is chosen
1014 INCL_DEBUG("OmegaN interaction: PiPiN channel chosen" << '\n');
1015 return new OmegaNToPiPiNChannel(particle1, particle2);
1016 }
1017 else {
1018 INCL_WARN("inconsistency within the OmegaN Cross Sections (sum!=inelastic)" << '\n');
1019 if(twoPiProductionCX>0.) {
1020 INCL_WARN("Returning a PiPiN channel" << '\n');
1021 isElastic = false;
1022 return new OmegaNToPiPiNChannel(particle1, particle2);
1023 } else if(onePiProductionCX>0.) {
1024 INCL_WARN("Returning a PiN channel" << '\n');
1025 isElastic = false;
1026 return new OmegaNToPiNChannel(particle1, particle2);
1027 } else {
1028 INCL_WARN("Returning an elastic channel" << '\n');
1029 isElastic = true;
1030 return new OmegaNElasticChannel(particle1, particle2);
1031 }
1032 }
1033 } else if ((particle1->isNucleon() && particle2->isKaon()) || (particle2->isNucleon() && particle1->isKaon())) {
1034//// KN
1035 const G4double elasticCX = CrossSections::elastic(particle1,particle2);
1036 const G4double quasielasticCX = CrossSections::NKToNK(particle1,particle2);
1037 const G4double NKToNKpiCX = CrossSections::NKToNKpi(particle1,particle2);
1038 const G4double NKToNK2piCX = CrossSections::NKToNK2pi(particle1,particle2);
1039 const G4double totCX=CrossSections::total(particle1, particle2);
1040// assert(std::fabs(totCX-elasticCX-quasielasticCX-NKToNKpiCX-NKToNK2piCX)<0.1);
1041
1042 const G4double rChannel=Random::shoot() * totCX;
1043 if(elasticCX > rChannel){
1044// Elastic KN channel is chosen
1045 isElastic = true;
1046 INCL_DEBUG("KN interaction: elastic channel chosen" << '\n');
1047 return new NKElasticChannel(particle1, particle2);
1048 } else if(elasticCX + quasielasticCX > rChannel){
1049// Quasi-elastic KN channel is chosen
1050 isElastic = false; // true ??
1051 INCL_DEBUG("KN interaction: quasi-elastic channel chosen" << '\n');
1052 return new NKToNKChannel(particle1, particle2);
1053 } else if(elasticCX + quasielasticCX + NKToNKpiCX > rChannel){
1054// KN -> NKpi channel is chosen
1055 isElastic = false;
1056 INCL_DEBUG("KN interaction: NKpi channel chosen" << '\n');
1057 return new NKToNKpiChannel(particle1, particle2);
1058 } else if(elasticCX + quasielasticCX + NKToNKpiCX + NKToNK2piCX > rChannel){
1059// KN -> NK2pi channel is chosen
1060 isElastic = false;
1061 INCL_DEBUG("KN interaction: NK2pi channel chosen" << '\n');
1062 return new NKToNK2piChannel(particle1, particle2);
1063 } else {
1064 INCL_WARN("inconsistency within the KN Cross Sections (sum!=inelastic)" << '\n');
1065 if(NKToNK2piCX>0.) {
1066 INCL_WARN("Returning a NKToNK2pi channel" << '\n');
1067 isElastic = false;
1068 return new NKToNK2piChannel(particle1, particle2);
1069 } else if(NKToNKpiCX>0.) {
1070 INCL_WARN("Returning a NKToNKpi channel" << '\n');
1071 isElastic = false;
1072 return new NKToNKpiChannel(particle1, particle2);
1073 } else if(quasielasticCX>0.) {
1074 INCL_WARN("Returning a quasi-elastic channel" << '\n');
1075 isElastic = false; // true ??
1076 return new NKToNKChannel(particle1, particle2);
1077 } else {
1078 INCL_WARN("Returning an elastic channel" << '\n');
1079 isElastic = true;
1080 return new NKElasticChannel(particle1, particle2);
1081 }
1082 }
1083 } else if ((particle1->isNucleon() && particle2->isAntiKaon()) || (particle2->isNucleon() && particle1->isAntiKaon())) {
1084//// KbN
1085 const G4double elasticCX = CrossSections::elastic(particle1,particle2);
1086 const G4double quasielasticCX = CrossSections::NKbToNKb(particle1,particle2);
1087 const G4double NKbToNKbpiCX = CrossSections::NKbToNKbpi(particle1,particle2);
1088 const G4double NKbToNKb2piCX = CrossSections::NKbToNKb2pi(particle1,particle2);
1089 const G4double NKbToLpiCX = CrossSections::NKbToLpi(particle1,particle2);
1090 const G4double NKbToL2piCX = CrossSections::NKbToL2pi(particle1,particle2);
1091 const G4double NKbToSpiCX = CrossSections::NKbToSpi(particle1,particle2);
1092 const G4double NKbToS2piCX = CrossSections::NKbToS2pi(particle1,particle2);
1093 const G4double totCX=CrossSections::total(particle1, particle2);
1094// assert(std::fabs(totCX-elasticCX-quasielasticCX-NKbToNKbpiCX-NKbToNKb2piCX-NKbToLpiCX-NKbToL2piCX-NKbToSpiCX-NKbToS2piCX)<0.1);
1095
1096 const G4double rChannel=Random::shoot() * totCX;
1097 if(elasticCX > rChannel){
1098// Elastic KbN channel is chosen
1099 isElastic = true;
1100 INCL_DEBUG("KbN interaction: elastic channel chosen" << '\n');
1101 return new NKbElasticChannel(particle1, particle2);
1102 } else if(elasticCX + quasielasticCX > rChannel){
1103// Quasi-elastic KbN channel is chosen
1104 isElastic = false; // true ??
1105 INCL_DEBUG("KbN interaction: quasi-elastic channel chosen" << '\n');
1106 return new NKbToNKbChannel(particle1, particle2);
1107 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX > rChannel){
1108// KbN -> NKbpi channel is chosen
1109 isElastic = false;
1110 INCL_DEBUG("KbN interaction: NKbpi channel chosen" << '\n');
1111 return new NKbToNKbpiChannel(particle1, particle2);
1112 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX + NKbToNKb2piCX > rChannel){
1113// KbN -> NKb2pi channel is chosen
1114 isElastic = false;
1115 INCL_DEBUG("KbN interaction: NKb2pi channel chosen" << '\n');
1116 return new NKbToNKb2piChannel(particle1, particle2);
1117 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX + NKbToNKb2piCX + NKbToLpiCX > rChannel){
1118// KbN -> Lpi channel is chosen
1119 isElastic = false;
1120 INCL_DEBUG("KbN interaction: Lpi channel chosen" << '\n');
1121 return new NKbToLpiChannel(particle1, particle2);
1122 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX + NKbToNKb2piCX + NKbToLpiCX + NKbToL2piCX > rChannel){
1123// KbN -> L2pi channel is chosen
1124 isElastic = false;
1125 INCL_DEBUG("KbN interaction: L2pi channel chosen" << '\n');
1126 return new NKbToL2piChannel(particle1, particle2);
1127 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX + NKbToNKb2piCX + NKbToLpiCX + NKbToL2piCX + NKbToSpiCX > rChannel){
1128// KbN -> Spi channel is chosen
1129 isElastic = false;
1130 INCL_DEBUG("KbN interaction: Spi channel chosen" << '\n');
1131 return new NKbToSpiChannel(particle1, particle2);
1132 } else if(elasticCX + quasielasticCX + NKbToNKbpiCX + NKbToNKb2piCX + NKbToLpiCX + NKbToL2piCX + NKbToSpiCX + NKbToS2piCX > rChannel){
1133// KbN -> S2pi channel is chosen
1134 isElastic = false;
1135 INCL_DEBUG("KbN interaction: S2pi channel chosen" << '\n');
1136 return new NKbToS2piChannel(particle1, particle2);
1137 } else {
1138 INCL_WARN("inconsistency within the KbN Cross Sections (sum!=inelastic)" << '\n');
1139 if(NKbToS2piCX>0.) {
1140 INCL_WARN("Returning a NKbToS2pi channel" << '\n');
1141 isElastic = false;
1142 return new NKbToS2piChannel(particle1, particle2);
1143 } else if(NKbToSpiCX>0.) {
1144 INCL_WARN("Returning a NKbToSpi channel" << '\n');
1145 isElastic = false;
1146 return new NKbToSpiChannel(particle1, particle2);
1147 } else if(NKbToL2piCX>0.) {
1148 INCL_WARN("Returning a NKbToL2pi channel" << '\n');
1149 isElastic = false;
1150 return new NKbToL2piChannel(particle1, particle2);
1151 } else if(NKbToLpiCX>0.) {
1152 INCL_WARN("Returning a NKbToLpi channel" << '\n');
1153 isElastic = false;
1154 return new NKbToLpiChannel(particle1, particle2);
1155 } else if(NKbToNKb2piCX>0.) {
1156 INCL_WARN("Returning a NKbToNKb2pi channel" << '\n');
1157 isElastic = false;
1158 return new NKbToNKb2piChannel(particle1, particle2);
1159 } else if(NKbToNKbpiCX>0.) {
1160 INCL_WARN("Returning a NKbToNKbpi channel" << '\n');
1161 isElastic = false;
1162 return new NKbToNKbpiChannel(particle1, particle2);
1163 } else if(quasielasticCX>0.) {
1164 INCL_WARN("Returning a quasi-elastic channel" << '\n');
1165 isElastic = false; // true ??
1166 return new NKbToNKbChannel(particle1, particle2);
1167 } else {
1168 INCL_WARN("Returning an elastic channel" << '\n');
1169 isElastic = true;
1170 return new NKbElasticChannel(particle1, particle2);
1171 }
1172 }
1173 } else if ((particle1->isNucleon() && particle2->isLambda()) || (particle2->isNucleon() && particle1->isLambda())) {
1174//// NLambda
1175 const G4double elasticCX = CrossSections::elastic(particle1,particle2);
1176 const G4double NLToNSCX = CrossSections::NLToNS(particle1,particle2);
1177 const G4double totCX=CrossSections::total(particle1, particle2);
1178// assert(std::fabs(totCX-elasticCX-NLToNSCX)<0.1);
1179
1180 const G4double rChannel=Random::shoot() * totCX;
1181 if(elasticCX > rChannel){
1182// Elastic NLambda channel is chosen
1183 isElastic = true;
1184 INCL_DEBUG("NLambda interaction: elastic channel chosen" << '\n');
1185 return new NYElasticChannel(particle1, particle2);
1186 } else if(elasticCX + NLToNSCX > rChannel){
1187// Quasi-elastic NLambda channel is chosen
1188 isElastic = false; // true ??
1189 INCL_DEBUG("NLambda interaction: quasi-elastic channel chosen" << '\n');
1190 return new NLToNSChannel(particle1, particle2);
1191 } else {
1192 INCL_WARN("inconsistency within the NLambda Cross Sections (sum!=inelastic)" << '\n');
1193 if(NLToNSCX>0.) {
1194 INCL_WARN("Returning a quasi-elastic channel" << '\n');
1195 isElastic = false; // true ??
1196 return new NLToNSChannel(particle1, particle2);
1197 } else {
1198 INCL_WARN("Returning an elastic channel" << '\n');
1199 isElastic = true;
1200 return new NYElasticChannel(particle1, particle2);
1201 }
1202 }
1203 } else if ((particle1->isNucleon() && particle2->isSigma()) || (particle2->isNucleon() && particle1->isSigma())) {
1204//// NSigma
1205 const G4double elasticCX = CrossSections::elastic(particle1,particle2);
1206 const G4double NSToNLCX = CrossSections::NSToNL(particle1,particle2);
1207 const G4double NSToNSCX = CrossSections::NSToNS(particle1,particle2);
1208 const G4double totCX=CrossSections::total(particle1, particle2);
1209// assert(std::fabs(totCX-elasticCX-NSToNLCX-NSToNSCX)<0.1);
1210
1211 const G4double rChannel=Random::shoot() * totCX;
1212 if(elasticCX > rChannel){
1213// Elastic NSigma channel is chosen
1214 isElastic = true;
1215 INCL_DEBUG("NSigma interaction: elastic channel chosen" << '\n');
1216 return new NYElasticChannel(particle1, particle2);
1217 } else if(elasticCX + NSToNLCX > rChannel){
1218// NSigma -> NLambda channel is chosen
1219 isElastic = false; // true ??
1220 INCL_DEBUG("NSigma interaction: NLambda channel chosen" << '\n');
1221 return new NSToNLChannel(particle1, particle2);
1222 } else if(elasticCX + NSToNLCX + NSToNSCX > rChannel){
1223// NSigma -> NSigma quasi-elastic channel is chosen
1224 isElastic = false; // true ??
1225 INCL_DEBUG("NSigma interaction: NSigma quasi-elastic channel chosen" << '\n');
1226 return new NSToNSChannel(particle1, particle2);
1227 } else {
1228 INCL_WARN("inconsistency within the NSigma Cross Sections (sum!=inelastic)" << '\n');
1229 if(NSToNSCX>0.) {
1230 INCL_WARN("Returning a quasi-elastic channel" << '\n');
1231 isElastic = false; // true ??
1232 return new NSToNSChannel(particle1, particle2);
1233 } else if(NSToNLCX>0.) {
1234 INCL_WARN("Returning a NLambda channel" << '\n');
1235 isElastic = false; // true ??
1236 return new NSToNLChannel(particle1, particle2);
1237 } else {
1238 INCL_WARN("Returning an elastic channel" << '\n');
1239 isElastic = true;
1240 return new NYElasticChannel(particle1, particle2);
1241 }
1242 }
1243 } else if ((particle1->isNucleon() && particle2->isAntiNucleon()) || (particle2->isNucleon() && particle1->isAntiNucleon())) {
1244//// NNbar
1245 const G4double totCX = CrossSections::total(particle1, particle2);
1246 const G4double NNbElasticCX = CrossSections::NNbarElastic(particle1,particle2);
1247 const G4double NNbCEXCX = CrossSections::NNbarCEX(particle1,particle2);
1248 const G4double NNbToLLbCX = CrossSections::NNbarToLLbar(particle1,particle2);
1249 const G4double NNbToNNbpiCX = CrossSections::NNbarToNNbarpi(particle1,particle2);
1250 const G4double NNbToNNb2piCX = CrossSections::NNbarToNNbar2pi(particle1,particle2);
1251 const G4double NNbToNNb3piCX = CrossSections::NNbarToNNbar3pi(particle1,particle2);
1252 const G4double AnnihilationCX = CrossSections::NNbarToAnnihilation(particle1, particle2);
1253
1254// assert(std::fabs(totCX-NNbElasticCX-NNbCEXCX-NNbToLLbCX-NNbToNNbpiCX-NNbToNNb2piCX-NNbToNNb3piCX-AnnihilationCX)<0.1);
1255
1256 const G4double rChannel=Random::shoot() * totCX;
1257 if (NNbElasticCX > rChannel) {
1258 // NNbar (elastic) channel is chosen
1259 isElastic = true;
1260 //INCL_WARN("NNbar interaction: NNbarElastic channel chosen" << '\n');
1261 return new NNbarElasticChannel(particle1, particle2);
1262 } else if (NNbElasticCX + NNbCEXCX > rChannel) {
1263 // NNbar (CEX) channel is chosen
1264 isElastic = false; // may be charge-exchange also
1265 //INCL_WARN("NNbar interaction: NNbarCEX channel chosen" << '\n');
1266 return new NNbarCEXChannel(particle1, particle2);
1267 } else if (NNbElasticCX + NNbCEXCX + NNbToLLbCX > rChannel) {
1268 // NNbarToLLbar channel is chosen
1269 isElastic = false; // may be charge-exchange also
1270 //INCL_WARN("NNbar interaction: NNbarToLLbar channel chosen" << '\n');
1271 return new NNbarToLLbarChannel(particle1, particle2);
1272 } else if (NNbElasticCX + NNbCEXCX + NNbToLLbCX + NNbToNNbpiCX > rChannel) {
1273 // NNbar to NNbar pi channel is chosen
1274 isElastic = false;
1275 //INCL_WARN("NNbar interaction: NNbar pi channel chosen" << '\n');
1276 return new NNbarToNNbarpiChannel(particle1, particle2);
1277 } else if (NNbElasticCX + NNbCEXCX + NNbToLLbCX + NNbToNNbpiCX + NNbToNNb2piCX > rChannel) {
1278 // NNbar to NNbar 2pi channel is chosen
1279 isElastic = false;
1280 //INCL_WARN("NNbar interaction: NNbar 2pi channel chosen" << '\n');
1281 return new NNbarToNNbar2piChannel(particle1, particle2);
1282 } else if (NNbElasticCX + NNbCEXCX + NNbToLLbCX + NNbToNNbpiCX + NNbToNNb2piCX + NNbToNNb3piCX > rChannel) {
1283 // NNbar to NNbar 3pi channel is chosen
1284 isElastic = false;
1285 //INCL_WARN("NNbar interaction: NNbar 3pi channel chosen" << '\n');
1286 return new NNbarToNNbar3piChannel(particle1, particle2);
1287 } else if (NNbElasticCX + NNbCEXCX + NNbToLLbCX + NNbToNNbpiCX + NNbToNNb2piCX + NNbToNNb3piCX +AnnihilationCX > rChannel){
1288 // NNbar annihilation channel is chosen
1289 isElastic = false;
1290 AnnihilationType atype;
1291 if((particle1->getType()==antiProton && particle2->getType()==Proton) || (particle2->getType()==antiProton && particle1->getType()==Proton)){
1292 atype = PTypeInFlight;
1293 }
1294 else if((particle1->getType()==antiProton && particle2->getType()==Neutron) || (particle2->getType()==antiProton && particle1->getType()==Neutron)){
1295 atype = NTypeInFlight;
1296 }
1297 else if((particle1->getType()==antiNeutron && particle2->getType()==Proton) || (particle2->getType()==antiNeutron && particle1->getType()==Proton)){
1298 atype = NbarPTypeInFlight;
1299 }
1300 else if((particle1->getType()==antiNeutron && particle2->getType()==Neutron) || (particle2->getType()==antiNeutron && particle1->getType()==Neutron)){
1301 atype = NbarNTypeInFlight;
1302 }
1303 else{
1304 atype = Def;
1305 INCL_ERROR("Annihilation type problem " << '\n');
1306 }
1307 theNucleus->setAType(atype);
1308 return new NNbarToAnnihilationChannel(theNucleus, particle1, particle2);
1309 } else {
1310 INCL_WARN("Inconsistency within the NNbar Cross Sections (sum != inelastic)" << '\n');
1311 if (NNbToNNb3piCX > 0.0) {
1312 INCL_WARN("Returning an NNbar 3pi channel" << '\n');
1313 isElastic = false;
1314 return new NNbarToNNbar3piChannel(particle1, particle2);
1315 } else if (NNbToNNb2piCX > 0.0) {
1316 INCL_WARN("Returning an NNbar 2pi channel" << '\n');
1317 isElastic = false;
1318 return new NNbarToNNbar2piChannel(particle1, particle2);
1319 } else if (NNbToNNbpiCX > 0.0) {
1320 INCL_WARN("Returning an NNbar pi channel" << '\n');
1321 isElastic = false;
1322 return new NNbarToNNbarpiChannel(particle1, particle2);
1323 } else if (AnnihilationCX > 0.0) {
1324 INCL_WARN("Returning an NNbar annihilation channel" << '\n');
1325 isElastic = false;
1326 AnnihilationType atype;
1327 if((particle1->getType()==antiProton && particle2->getType()==Proton) || (particle2->getType()==antiProton && particle1->getType()==Proton)){
1328 atype = PTypeInFlight;
1329 }
1330 else if((particle1->getType()==antiProton && particle2->getType()==Neutron) || (particle2->getType()==antiProton && particle1->getType()==Neutron)){
1331 atype = NTypeInFlight;
1332 }
1333 else if((particle1->getType()==antiNeutron && particle2->getType()==Proton) || (particle2->getType()==antiNeutron && particle1->getType()==Proton)){
1334 atype = NbarPTypeInFlight;
1335 }
1336 else if((particle1->getType()==antiNeutron && particle2->getType()==Neutron) || (particle2->getType()==antiNeutron && particle1->getType()==Neutron)){
1337 atype = NbarNTypeInFlight;
1338 }
1339 else{
1340 atype = Def;
1341 INCL_ERROR("Annihilation type problem " << '\n');
1342 }
1343 theNucleus->setAType(atype);
1344 return new NNbarToAnnihilationChannel(theNucleus, particle1, particle2);
1345 } else if (NNbCEXCX > 0.0) {
1346 INCL_WARN("Returning an NNbar CEX channel" << '\n');
1347 isElastic = false;
1348 return new NNbarCEXChannel(particle1, particle2);
1349 } else if (NNbToLLbCX > 0.0) {
1350 INCL_WARN("Returning an NNbar LLbar channel" << '\n');
1351 isElastic = false;
1352 return new NNbarToLLbarChannel(particle1, particle2);
1353 } else {
1354 INCL_WARN("Elastic NNbar channel chosen" << '\n');
1355 isElastic = true;
1356 return new NNbarElasticChannel(particle1, particle2);
1357 }
1358 }
1359 }
1360
1361 else {
1362 INCL_DEBUG("BinaryCollisionAvatar can only handle nucleons (for the moment)."
1363 << '\n'
1364 << particle1->print()
1365 << '\n'
1366 << particle2->print()
1367 << '\n');
1368 InteractionAvatar::restoreParticles();
1369 return NULL;
1370 }
1371 }
INCL_ERROR
#define INCL_ERROR(x)
Definition G4INCLLogger.hh:225
INCL_WARN
#define INCL_WARN(x)
Definition G4INCLLogger.hh:232
INCL_DEBUG
#define INCL_DEBUG(x)
Definition G4INCLLogger.hh:240
G4double
double G4double
Definition G4Types.hh:83
G4INCL::Book::getAcceptedCollisions
G4int getAcceptedCollisions() const
Definition G4INCLBook.hh:100
G4INCL::InteractionAvatar::restoreParticles
void restoreParticles() const
Restore the state of both particles.
Definition G4INCLInteractionAvatar.cc:341
G4INCL::Nucleus::getStore
Store * getStore() const
Definition G4INCLNucleus.hh:361
G4INCL::Nucleus::setAType
void setAType(AnnihilationType type)
Definition G4INCLNucleus.cc:133
G4INCL::Particle::isLambda
G4bool isLambda() const
Is this a Lambda?
Definition G4INCLParticle.hh:434
G4INCL::Particle::isOmega
G4bool isOmega() const
Is this an omega?
Definition G4INCLParticle.hh:408
G4INCL::Particle::isSigma
G4bool isSigma() const
Is this a Sigma?
Definition G4INCLParticle.hh:425
G4INCL::Particle::getPosition
const G4INCL::ThreeVector & getPosition() const
Definition G4INCLParticle.hh:956
G4INCL::Particle::MergeVectorBias
static std::vector< G4int > MergeVectorBias(Particle const *const p1, Particle const *const p2)
Definition G4INCLParticle.cc:235
G4INCL::Particle::isAntiNucleon
G4bool isAntiNucleon() const
Is this an antinucleon?
Definition G4INCLParticle.hh:455
G4INCL::Particle::isEta
G4bool isEta() const
Is this an eta?
Definition G4INCLParticle.hh:405
G4INCL::Particle::isAntiKaon
G4bool isAntiKaon() const
Is this an antiKaon?
Definition G4INCLParticle.hh:431
G4INCL::Particle::getType
G4INCL::ParticleType getType() const
Definition G4INCLParticle.hh:187
G4INCL::Particle::isKaon
G4bool isKaon() const
Is this a Kaon?
Definition G4INCLParticle.hh:428
G4INCL::Particle::print
std::string print() const
Definition G4INCLParticle.hh:1093
G4INCL::Particle::isDelta
G4bool isDelta() const
Is it a Delta?
Definition G4INCLParticle.hh:420
G4INCL::Particle::getBiasFromVector
static G4double getBiasFromVector(std::vector< G4int > VectorBias)
Definition G4INCLParticle.cc:223
G4INCL::Particle::isNucleon
G4bool isNucleon() const
Definition G4INCLParticle.hh:362
G4INCL::Store::getBook
Book & getBook()
Definition G4INCLStore.hh:259
G4INCL::ThreeVector::dot
G4double dot(const ThreeVector &v) const
Definition G4INCLThreeVector.hh:97
G4INCL::ThreeVector::mag2
G4double mag2() const
Definition G4INCLThreeVector.hh:78
G4INCL::CrossSections::NSToNL
G4double NSToNL(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:260
G4INCL::CrossSections::NNToNNKKb
G4double NNToNNKKb(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:187
G4INCL::CrossSections::elastic
G4double elastic(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:57
G4INCL::CrossSections::NpiToSK
G4double NpiToSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:216
G4INCL::CrossSections::piNToOmegaN
G4double piNToOmegaN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:89
G4INCL::CrossSections::NNbarCEX
G4double NNbarCEX(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:295
G4INCL::CrossSections::NKbToNKb2pi
G4double NKbToNKb2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:328
G4INCL::CrossSections::NNToNDeltaOmega
G4double NNToNDeltaOmega(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:146
G4INCL::CrossSections::NDeltaToDeltaLK
G4double NDeltaToDeltaLK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:201
G4INCL::CrossSections::NNToNSKpi
G4double NNToNSKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:175
G4INCL::CrossSections::etaNToPiN
G4double etaNToPiN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:97
G4INCL::CrossSections::NDeltaToNNKKb
G4double NDeltaToNNKKb(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:208
G4INCL::CrossSections::NKbToLpi
G4double NKbToLpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:288
G4INCL::CrossSections::NNToNLK2pi
G4double NNToNLK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:179
G4INCL::CrossSections::etaNToPiPiN
G4double etaNToPiPiN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:101
G4INCL::CrossSections::NKbToS2pi
G4double NKbToS2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:316
G4INCL::CrossSections::piNToEtaN
G4double piNToEtaN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:85
G4INCL::CrossSections::NNbarToAnnihilation
G4double NNbarToAnnihilation(Particle const *const p1, Particle const *const p2)
Nucleon-AntiNucleon total annihilation cross sections.
Definition G4INCLCrossSections.cc:312
G4INCL::CrossSections::omegaNToPiN
G4double omegaNToPiN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:105
G4INCL::CrossSections::NDeltaToNSK
G4double NDeltaToNSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:198
G4INCL::CrossSections::NKbToSpi
G4double NKbToSpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:284
G4INCL::CrossSections::NNToNDelta
G4double NNToNDelta(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:69
G4INCL::CrossSections::NNToNSK
G4double NNToNSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:167
G4INCL::CrossSections::NDeltaToNLK
G4double NDeltaToNLK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:195
G4INCL::CrossSections::NLToNS
G4double NLToNS(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:256
G4INCL::CrossSections::piNToDelta
G4double piNToDelta(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:77
G4INCL::CrossSections::NKbToNKb
G4double NKbToNKb(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:280
G4INCL::CrossSections::NNToNSK2pi
G4double NNToNSK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:183
G4INCL::CrossSections::NNbarToNNbar3pi
G4double NNbarToNNbar3pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:308
G4INCL::CrossSections::NNToNLK
G4double NNToNLK(Particle const *const p1, Particle const *const p2)
Strange cross sections.
Definition G4INCLCrossSections.cc:163
G4INCL::CrossSections::NNToNNEtaExclu
G4double NNToNNEtaExclu(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:121
G4INCL::CrossSections::NNToNDeltaEta
G4double NNToNDeltaEta(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:129
G4INCL::CrossSections::NNToNNOmegaxPi
G4double NNToNNOmegaxPi(const G4int xpi, Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:142
G4INCL::CrossSections::NNToNNEtaxPi
G4double NNToNNEtaxPi(const G4int xpi, Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:125
G4INCL::CrossSections::NNbarToNNbarpi
G4double NNbarToNNbarpi(Particle const *const p1, Particle const *const p2)
Nucleon-AntiNucleon to Nucleon-AntiNucleon + pions cross sections.
Definition G4INCLCrossSections.cc:302
G4INCL::CrossSections::NSToNS
G4double NSToNS(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:264
G4INCL::CrossSections::NpiToNKKb
G4double NpiToNKKb(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:248
G4INCL::CrossSections::NpiToLK2pi
G4double NpiToLK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:240
G4INCL::CrossSections::omegaNToPiPiN
G4double omegaNToPiPiN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:109
G4INCL::CrossSections::NDeltaToDeltaSK
G4double NDeltaToDeltaSK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:204
G4INCL::CrossSections::NpiToMissingStrangeness
G4double NpiToMissingStrangeness(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:252
G4INCL::CrossSections::NNbarToNNbar2pi
G4double NNbarToNNbar2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:305
G4INCL::CrossSections::total
G4double total(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:61
G4INCL::CrossSections::NDeltaToNN
G4double NDeltaToNN(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:65
G4INCL::CrossSections::NpiToLKpi
G4double NpiToLKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:232
G4INCL::CrossSections::NKbToNKbpi
G4double NKbToNKbpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:324
G4INCL::CrossSections::NNToNLKpi
G4double NNToNLKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:171
G4INCL::CrossSections::NKbToL2pi
G4double NKbToL2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:320
G4INCL::CrossSections::NNToxPiNN
G4double NNToxPiNN(const G4int xpi, Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:73
G4INCL::CrossSections::NKToNK2pi
G4double NKToNK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:276
G4INCL::CrossSections::NKToNKpi
G4double NKToNKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:272
G4INCL::CrossSections::NNToNNOmegaExclu
G4double NNToNNOmegaExclu(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:138
G4INCL::CrossSections::NpiToSK2pi
G4double NpiToSK2pi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:244
G4INCL::CrossSections::NKToNK
G4double NKToNK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:268
G4INCL::CrossSections::NNToMissingStrangeness
G4double NNToMissingStrangeness(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:191
G4INCL::CrossSections::piNToxPiN
G4double piNToxPiN(const G4int xpi, Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:81
G4INCL::CrossSections::NNbarToLLbar
G4double NNbarToLLbar(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:298
G4INCL::CrossSections::NpiToLK
G4double NpiToLK(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:212
G4INCL::CrossSections::NNbarElastic
G4double NNbarElastic(Particle const *const p1, Particle const *const p2)
antiparticle cross sections
Definition G4INCLCrossSections.cc:292
G4INCL::CrossSections::NpiToSKpi
G4double NpiToSKpi(Particle const *const p1, Particle const *const p2)
Definition G4INCLCrossSections.cc:236
G4INCL::KinematicsUtils::squareTotalEnergyInCM
G4double squareTotalEnergyInCM(Particle const *const p1, Particle const *const p2)
Definition G4INCLKinematicsUtils.cc:125
G4INCL::Math::tenPi
const G4double tenPi
Definition G4INCLGlobals.hh:70
G4INCL::Random::shoot
G4double shoot()
Definition G4INCLRandom.cc:93
G4INCL::NbarPTypeInFlight
@ NbarPTypeInFlight
Definition G4INCLNucleus.hh:65
G4INCL::NbarNTypeInFlight
@ NbarNTypeInFlight
Definition G4INCLNucleus.hh:65

◆ getCutNN()

static G4double G4INCL::BinaryCollisionAvatar::getCutNN ( )
inlinestatic

Definition at line 79 of file G4INCLBinaryCollisionAvatar.hh.

79{ return cutNN; }

◆ getCutNNSquared()

static G4double G4INCL::BinaryCollisionAvatar::getCutNNSquared ( )
inlinestatic

Definition at line 81 of file G4INCLBinaryCollisionAvatar.hh.

81{ return cutNNSquared; }

Referenced by G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar().

◆ getParticles()

ParticleList G4INCL::BinaryCollisionAvatar::getParticles ( ) const
inlinevirtual

Implements G4INCL::IAvatar.

Definition at line 62 of file G4INCLBinaryCollisionAvatar.hh.

62 {
63 ParticleList theParticleList;
64 theParticleList.push_back(particle1);
65 theParticleList.push_back(particle2);
66 return theParticleList;
67 };

◆ postInteraction()

void G4INCL::BinaryCollisionAvatar::postInteraction ( FinalState * fs)
virtual

Implements G4INCL::IAvatar.

Definition at line 1379 of file G4INCLBinaryCollisionAvatar.cc.

1379 {
1380 // Call the postInteraction method of the parent class
1381 // (provides Pauli blocking and enforces energy conservation)
1382 InteractionAvatar::postInteraction(fs);
1383
1384 switch(fs->getValidity()) {
1385 case PauliBlockedFS:
1386 theNucleus->getStore()->getBook().incrementBlockedCollisions();
1387 break;
1388 case NoEnergyConservationFS:
1389 case ParticleBelowFermiFS:
1390 case ParticleBelowZeroFS:
1391 break;
1392 case ValidFS:
1393 Book &theBook = theNucleus->getStore()->getBook();
1394 theBook.incrementAcceptedCollisions();
1395
1396 if(theBook.getAcceptedCollisions() == 1) {
1397 // Store time and cross section of the first collision
1398 G4double t = theBook.getCurrentTime();
1399 theBook.setFirstCollisionTime(t);
1400 theBook.setFirstCollisionXSec(oldXSec);
1401 // Increase the number of Kaon by 1
1402 if(isStrangeProduction) theNucleus->setNumberOfKaon(theNucleus->getNumberOfKaon()+1);
1403 // Store position and momentum of the spectator on the first
1404 // collision
1405 if((isParticle1Spectator && isParticle2Spectator) || (!isParticle1Spectator && !isParticle2Spectator)) {
1406 INCL_ERROR("First collision must be within a target spectator and a non-target spectator");
1407 }
1408 if(isParticle1Spectator) {
1409 theBook.setFirstCollisionSpectatorPosition(backupParticle1->getPosition().mag());
1410 theBook.setFirstCollisionSpectatorMomentum(backupParticle1->getMomentum().mag());
1411 } else {
1412 theBook.setFirstCollisionSpectatorPosition(backupParticle2->getPosition().mag());
1413 theBook.setFirstCollisionSpectatorMomentum(backupParticle2->getMomentum().mag());
1414 }
1415
1416 // Store the elasticity of the first collision
1417 theBook.setFirstCollisionIsElastic(isElastic);
1418 }
1419 }
1420 return;
1421 }
G4INCL::Book::incrementAcceptedCollisions
void incrementAcceptedCollisions()
Definition G4INCLBook.hh:72
G4INCL::Book::incrementBlockedCollisions
void incrementBlockedCollisions()
Definition G4INCLBook.hh:73
G4INCL::InteractionAvatar::backupParticle2
static G4ThreadLocal Particle * backupParticle2
Definition G4INCLInteractionAvatar.hh:108
G4INCL::InteractionAvatar::backupParticle1
static G4ThreadLocal Particle * backupParticle1
Definition G4INCLInteractionAvatar.hh:108
G4INCL::Particle::setNumberOfKaon
void setNumberOfKaon(const G4int NK)
Definition G4INCLParticle.hh:1194
G4INCL::Particle::getMomentum
const G4INCL::ThreeVector & getMomentum() const
Definition G4INCLParticle.hh:934
G4INCL::Particle::getNumberOfKaon
G4int getNumberOfKaon() const
Number of Kaon inside de nucleus.
Definition G4INCLParticle.hh:1193
G4INCL::NoEnergyConservationFS
@ NoEnergyConservationFS
Definition G4INCLFinalState.hh:56

◆ preInteraction()

void G4INCL::BinaryCollisionAvatar::preInteraction ( )
virtual

Implements G4INCL::IAvatar.

Definition at line 1373 of file G4INCLBinaryCollisionAvatar.cc.

1373 {
1374 isParticle1Spectator = particle1->isTargetSpectator();
1375 isParticle2Spectator = particle2->isTargetSpectator();
1376 InteractionAvatar::preInteraction();
1377 }
G4INCL::Particle::isTargetSpectator
G4bool isTargetSpectator() const
Definition G4INCLParticle.hh:381

◆ setBias()

static void G4INCL::BinaryCollisionAvatar::setBias ( const G4double b)
inlinestatic

Set the global bias factor.

Definition at line 87 of file G4INCLBinaryCollisionAvatar.hh.

87{ bias=b; }

Referenced by G4INCL::INCL::INCL().

◆ setCutNN()

static void G4INCL::BinaryCollisionAvatar::setCutNN ( const G4double c)
inlinestatic

Definition at line 74 of file G4INCLBinaryCollisionAvatar.hh.

74 {
75 cutNN = c;
76 cutNNSquared = cutNN*cutNN;
77 }

Referenced by G4INCL::INCL::INCL().

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