Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4WilsonAbrasionModel Class Reference

#include <G4WilsonAbrasionModel.hh>

Inheritance diagram for G4WilsonAbrasionModel:
Collaboration diagram for G4WilsonAbrasionModel:

Public Member Functions

 G4WilsonAbrasionModel (G4bool useAblation1=false)
 
 G4WilsonAbrasionModel (G4ExcitationHandler *)
 
 ~G4WilsonAbrasionModel ()
 
 G4WilsonAbrasionModel (const G4WilsonAbrasionModel &right)
 
const G4WilsonAbrasionModeloperator= (G4WilsonAbrasionModel &right)
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &, G4Nucleus &)
 
void SetVerboseLevel (G4int)
 
void SetUseAblation (G4bool)
 
G4bool GetUseAblation ()
 
void SetConserveMomentum (G4bool)
 
G4bool GetConserveMomentum ()
 
void SetExcitationHandler (G4ExcitationHandler *)
 
G4ExcitationHandlerGetExcitationHandler ()
 
virtual void ModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double
GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double
GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 77 of file G4WilsonAbrasionModel.hh.

Constructor & Destructor Documentation

G4WilsonAbrasionModel::G4WilsonAbrasionModel ( G4bool  useAblation1 = false)

Definition at line 117 of file G4WilsonAbrasionModel.cc.

118  :G4HadronicInteraction("G4WilsonAbrasion")
119 {
120  // Send message to stdout to advise that the G4Abrasion model is being used.
121  PrintWelcomeMessage();
122 
123  // Set the default verbose level to 0 - no output.
124  verboseLevel = 0;
125  useAblation = useAblation1;
126 
127  // No de-excitation handler has been supplied - define the default handler.
128 
129  theExcitationHandler = new G4ExcitationHandler;
130  theExcitationHandlerx = new G4ExcitationHandler;
131  if (useAblation)
132  {
133  theAblation = new G4WilsonAblationModel;
134  theAblation->SetVerboseLevel(verboseLevel);
135  theExcitationHandler->SetEvaporation(theAblation);
136  theExcitationHandlerx->SetEvaporation(theAblation);
137  }
138  else
139  {
140  theAblation = NULL;
141  G4Evaporation * theEvaporation = new G4Evaporation;
142  G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
143  G4StatMF * theMF = new G4StatMF;
144  theExcitationHandler->SetEvaporation(theEvaporation);
145  theExcitationHandler->SetFermiModel(theFermiBreakUp);
146  theExcitationHandler->SetMultiFragmentation(theMF);
147  theExcitationHandler->SetMaxAandZForFermiBreakUp(12, 6);
148  theExcitationHandler->SetMinEForMultiFrag(5.0*MeV);
149 
150  theEvaporation = new G4Evaporation;
151  theFermiBreakUp = new G4FermiBreakUp;
152  theExcitationHandlerx->SetEvaporation(theEvaporation);
153  theExcitationHandlerx->SetFermiModel(theFermiBreakUp);
154  theExcitationHandlerx->SetMaxAandZForFermiBreakUp(12, 6);
155  }
156 
157  // Set the minimum and maximum range for the model (despite nomanclature,
158  // this is in energy per nucleon number).
159 
160  SetMinEnergy(70.0*MeV);
161  SetMaxEnergy(10.1*GeV);
162  isBlocked = false;
163 
164  // npK, when mutiplied by the nuclear Fermi momentum, determines the range of
165  // momentum over which the secondary nucleon momentum is sampled.
166 
167  r0sq = 0.0;
168  npK = 5.0;
169  B = 10.0 * MeV;
170  third = 1.0 / 3.0;
171  fradius = 0.99;
172  conserveEnergy = false;
173  conserveMomentum = true;
174 }
void SetMinEForMultiFrag(G4double anE)
void SetMinEnergy(G4double anEnergy)
void SetFermiModel(G4VFermiBreakUp *ptr)
void SetMultiFragmentation(G4VMultiFragmentation *ptr)
G4HadronicInteraction(const G4String &modelName="HadronicModel")
void SetEvaporation(G4VEvaporation *ptr, G4bool isLocal=false)
void SetMaxAandZForFermiBreakUp(G4int anA, G4int aZ)
static constexpr double GeV
Definition: G4SIunits.hh:217
void SetMaxEnergy(const G4double anEnergy)
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

G4WilsonAbrasionModel::G4WilsonAbrasionModel ( G4ExcitationHandler aExcitationHandler)

Definition at line 188 of file G4WilsonAbrasionModel.cc.

189 {
190 // Send message to stdout to advise that the G4Abrasion model is being used.
191 
192  PrintWelcomeMessage();
193 
194 // Set the default verbose level to 0 - no output.
195 
196  verboseLevel = 0;
197 
198  theAblation = NULL; //A.R. 26-Jul-2012 Coverity fix.
199  useAblation = false; //A.R. 14-Aug-2012 Coverity fix.
200 
201 //
202 // The user is able to provide the excitation handler as well as an argument
203 // which is provided in this instantiation is used to determine
204 // whether the spectators of the interaction are free following the abrasion.
205 //
206  theExcitationHandler = aExcitationHandler;
207  theExcitationHandlerx = new G4ExcitationHandler;
208  G4Evaporation * theEvaporation = new G4Evaporation;
209  G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
210  theExcitationHandlerx->SetEvaporation(theEvaporation);
211  theExcitationHandlerx->SetFermiModel(theFermiBreakUp);
212  theExcitationHandlerx->SetMaxAandZForFermiBreakUp(12, 6);
213 //
214 //
215 // Set the minimum and maximum range for the model (despite nomanclature, this
216 // is in energy per nucleon number).
217 //
218  SetMinEnergy(70.0*MeV);
219  SetMaxEnergy(10.1*GeV);
220  isBlocked = false;
221 //
222 //
223 // npK, when mutiplied by the nuclear Fermi momentum, determines the range of
224 // momentum over which the secondary nucleon momentum is sampled.
225 //
226  r0sq = 0.0; //A.R. 14-Aug-2012 Coverity fix.
227  npK = 5.0;
228  B = 10.0 * MeV;
229  third = 1.0 / 3.0;
230  fradius = 0.99;
231  conserveEnergy = false;
232  conserveMomentum = true;
233 }
void SetMinEnergy(G4double anEnergy)
void SetFermiModel(G4VFermiBreakUp *ptr)
void SetEvaporation(G4VEvaporation *ptr, G4bool isLocal=false)
void SetMaxAandZForFermiBreakUp(G4int anA, G4int aZ)
static constexpr double GeV
Definition: G4SIunits.hh:217
void SetMaxEnergy(const G4double anEnergy)
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

G4WilsonAbrasionModel::~G4WilsonAbrasionModel ( )

Definition at line 236 of file G4WilsonAbrasionModel.cc.

237 {
238 //
239 //
240 // The destructor doesn't have to do a great deal!
241 //
242  if (theExcitationHandler) delete theExcitationHandler;
243  if (theExcitationHandlerx) delete theExcitationHandlerx;
244  if (useAblation) delete theAblation;
245 // delete theExcitationHandler;
246 // delete theExcitationHandlerx;
247 }
G4WilsonAbrasionModel::G4WilsonAbrasionModel ( const G4WilsonAbrasionModel right)

Member Function Documentation

G4HadFinalState * G4WilsonAbrasionModel::ApplyYourself ( const G4HadProjectile theTrack,
G4Nucleus theTarget 
)
virtual

Implements G4HadronicInteraction.

Definition at line 250 of file G4WilsonAbrasionModel.cc.

252 {
253 //
254 //
255 // The secondaries will be returned in G4HadFinalState &theParticleChange -
256 // initialise this. The original track will always be discontinued and
257 // secondaries followed.
258 //
261 //
262 //
263 // Get relevant information about the projectile and target (A, Z, energy/nuc,
264 // momentum, etc).
265 //
266  const G4ParticleDefinition *definitionP = theTrack.GetDefinition();
267  const G4double AP = definitionP->GetBaryonNumber();
268  const G4double ZP = definitionP->GetPDGCharge();
269  G4LorentzVector pP = theTrack.Get4Momentum();
270  G4double E = theTrack.GetKineticEnergy()/AP;
271  G4double AT = theTarget.GetA_asInt();
272  G4double ZT = theTarget.GetZ_asInt();
273  G4double TotalEPre = theTrack.GetTotalEnergy() +
274  theTarget.AtomicMass(AT, ZT) + theTarget.GetEnergyDeposit();
275  G4double TotalEPost = 0.0;
276 //
277 //
278 // Determine the radii of the projectile and target nuclei.
279 //
280  G4WilsonRadius aR;
281  G4double rP = aR.GetWilsonRadius(AP);
282  G4double rT = aR.GetWilsonRadius(AT);
283  G4double rPsq = rP * rP;
284  G4double rTsq = rT * rT;
285  if (verboseLevel >= 2)
286  {
287  G4cout <<"########################################"
288  <<"########################################"
289  <<G4endl;
290  G4cout.precision(6);
291  G4cout <<"IN G4WilsonAbrasionModel" <<G4endl;
292  G4cout <<"Initial projectile A=" <<AP
293  <<", Z=" <<ZP
294  <<", radius = " <<rP/fermi <<" fm"
295  <<G4endl;
296  G4cout <<"Initial target A=" <<AT
297  <<", Z=" <<ZT
298  <<", radius = " <<rT/fermi <<" fm"
299  <<G4endl;
300  G4cout <<"Projectile momentum and Energy/nuc = " <<pP <<" ," <<E <<G4endl;
301  }
302 //
303 //
304 // The following variables are used to determine the impact parameter in the
305 // near-field (i.e. taking into consideration the electrostatic repulsion).
306 //
307  G4double rm = ZP * ZT * elm_coupling / (E * AP);
308  G4double r = 0.0;
309  G4double rsq = 0.0;
310 //
311 //
312 // Initialise some of the variables which wll be used to calculate the chord-
313 // length for nucleons in the projectile and target, and hence calculate the
314 // number of abraded nucleons and the excitation energy.
315 //
316  G4NuclearAbrasionGeometry *theAbrasionGeometry = NULL;
317  G4double CT = 0.0;
318  G4double F = 0.0;
319  G4int Dabr = 0;
320 //
321 //
322 // The following loop is performed until the number of nucleons which are
323 // abraded by the process is >1, i.e. an interaction MUST occur.
324 //
325  G4bool skipInteraction = false; // It will be set true if the two nuclei fail to collide
326  const G4int maxNumberOfLoops = 1000;
327  G4int loopCounter = -1;
328  while (Dabr == 0 && ++loopCounter < maxNumberOfLoops) /* Loop checking, 07.08.2015, A.Ribon */
329  {
330 // Added by MHM 20050119 to fix leaking memory on second pass through this loop
331  if (theAbrasionGeometry)
332  {
333  delete theAbrasionGeometry;
334  theAbrasionGeometry = NULL;
335  }
336 //
337 //
338 // Sample the impact parameter. For the moment, this class takes account of
339 // electrostatic effects on the impact parameter, but (like HZETRN AND NUCFRG2)
340 // does not make any correction for the effects of nuclear-nuclear repulsion.
341 //
342  G4double rPT = rP + rT;
343  G4double rPTsq = rPT * rPT;
344 //
345 //
346 // This is a "catch" to make sure we don't go into an infinite loop because the
347 // energy is too low to overcome nuclear repulsion. PRT 20091023. If the
348 // value of rm < fradius * rPT then we're unlikely to sample a small enough
349 // impact parameter (energy of incident particle is too low).
350 //
351  if (rm >= fradius * rPT) {
352  skipInteraction = true;
353  }
354 //
355 //
356 // Now sample impact parameter until the criterion is met that projectile
357 // and target overlap, but repulsion is taken into consideration.
358 //
359  G4int evtcnt = 0;
360  r = 1.1 * rPT;
361  while (r > rPT && ++evtcnt < 1000) /* Loop checking, 07.08.2015, A.Ribon */
362  {
363  G4double bsq = rPTsq * G4UniformRand();
364  r = (rm + std::sqrt(rm*rm + 4.0*bsq)) / 2.0;
365  }
366 //
367 //
368 // We've tried to sample this 1000 times, but failed.
369 //
370  if (evtcnt >= 1000) {
371  skipInteraction = true;
372  }
373 
374  rsq = r * r;
375 //
376 //
377 // Now determine the chord-length through the target nucleus.
378 //
379  if (rT > rP)
380  {
381  G4double x = (rPsq + rsq - rTsq) / 2.0 / r;
382  if (x > 0.0) CT = 2.0 * std::sqrt(rTsq - x*x);
383  else CT = 2.0 * std::sqrt(rTsq - rsq);
384  }
385  else
386  {
387  G4double x = (rTsq + rsq - rPsq) / 2.0 / r;
388  if (x > 0.0) CT = 2.0 * std::sqrt(rTsq - x*x);
389  else CT = 2.0 * rT;
390  }
391 //
392 //
393 // Determine the number of abraded nucleons. Note that the mean number of
394 // abraded nucleons is used to sample the Poisson distribution. The Poisson
395 // distribution is sampled only ten times with the current impact parameter,
396 // and if it fails after this to find a case for which the number of abraded
397 // nucleons >1, the impact parameter is re-sampled.
398 //
399  theAbrasionGeometry = new G4NuclearAbrasionGeometry(AP,AT,r);
400  F = theAbrasionGeometry->F();
401  G4double lambda = 16.6*fermi / G4Pow::GetInstance()->powA(E/MeV,0.26);
402  G4double Mabr = F * AP * (1.0 - G4Exp(-CT/lambda));
403  G4long n = 0;
404  for (G4int i = 0; i<10; i++)
405  {
406  n = G4Poisson(Mabr);
407  if (n > 0)
408  {
409  if (n>AP) Dabr = (G4int) AP;
410  else Dabr = (G4int) n;
411  break;
412  }
413  }
414  } // End of while loop
415 
416  if ( loopCounter >= maxNumberOfLoops || skipInteraction ) {
417  // Assume nuclei do not collide and return unchanged primary.
421  if (verboseLevel >= 2) {
422  G4cout <<"Particle energy too low to overcome repulsion." <<G4endl;
423  G4cout <<"Event rejected and original track maintained" <<G4endl;
424  G4cout <<"########################################"
425  <<"########################################"
426  <<G4endl;
427  }
428  return &theParticleChange;
429  }
430 
431  if (verboseLevel >= 2)
432  {
433  G4cout <<G4endl;
434  G4cout <<"Impact parameter = " <<r/fermi <<" fm" <<G4endl;
435  G4cout <<"# Abraded nucleons = " <<Dabr <<G4endl;
436  }
437 //
438 //
439 // The number of abraded nucleons must be no greater than the number of
440 // nucleons in either the projectile or the target. If AP - Dabr < 2 or
441 // AT - Dabr < 2 then either we have only a nucleon left behind in the
442 // projectile/target or we've tried to abrade too many nucleons - and Dabr
443 // should be limited.
444 //
445  if (AP - (G4double) Dabr < 2.0) Dabr = (G4int) AP;
446  if (AT - (G4double) Dabr < 2.0) Dabr = (G4int) AT;
447 //
448 //
449 // Determine the abraded secondary nucleons from the projectile. *fragmentP
450 // is a pointer to the prefragment from the projectile and nSecP is the number
451 // of nucleons in theParticleChange which have been abraded. The total energy
452 // from these is determined.
453 //
454  G4ThreeVector boost = pP.findBoostToCM();
455  G4Fragment *fragmentP = GetAbradedNucleons (Dabr, AP, ZP, rP);
457  G4int i = 0;
458  for (i=0; i<nSecP; i++)
459  {
460  TotalEPost += theParticleChange.GetSecondary(i)->
461  GetParticle()->GetTotalEnergy();
462  }
463 //
464 //
465 // Determine the number of spectators in the interaction region for the
466 // projectile.
467 //
468  G4int DspcP = (G4int) (AP*F) - Dabr;
469  if (DspcP <= 0) DspcP = 0;
470  else if (DspcP > AP-Dabr) DspcP = ((G4int) AP) - Dabr;
471 //
472 //
473 // Determine excitation energy associated with excess surface area of the
474 // projectile (EsP) and the excitation due to scattering of nucleons which are
475 // retained within the projectile (ExP). Add the total energy from the excited
476 // nucleus to the total energy of the secondaries.
477 //
478  G4bool excitationAbsorbedByProjectile = false;
479  if (fragmentP != NULL)
480  {
481  G4double EsP = theAbrasionGeometry->GetExcitationEnergyOfProjectile();
482  G4double ExP = 0.0;
483  if (Dabr < AT)
484  excitationAbsorbedByProjectile = G4UniformRand() < 0.5;
485  if (excitationAbsorbedByProjectile)
486  ExP = GetNucleonInducedExcitation(rP, rT, r);
487  G4double xP = EsP + ExP;
488  if (xP > B*(AP-Dabr)) xP = B*(AP-Dabr);
489  G4LorentzVector lorentzVector = fragmentP->GetMomentum();
490  lorentzVector.setE(lorentzVector.e()+xP);
491  fragmentP->SetMomentum(lorentzVector);
492  TotalEPost += lorentzVector.e();
493  }
494  G4double EMassP = TotalEPost;
495 //
496 //
497 // Determine the abraded secondary nucleons from the target. Note that it's
498 // assumed that the same number of nucleons are abraded from the target as for
499 // the projectile, and obviously no boost is applied to the products. *fragmentT
500 // is a pointer to the prefragment from the target and nSec is the total number
501 // of nucleons in theParticleChange which have been abraded. The total energy
502 // from these is determined.
503 //
504  G4Fragment *fragmentT = GetAbradedNucleons (Dabr, AT, ZT, rT);
506  for (i=nSecP; i<nSec; i++)
507  {
508  TotalEPost += theParticleChange.GetSecondary(i)->
509  GetParticle()->GetTotalEnergy();
510  }
511 //
512 //
513 // Determine the number of spectators in the interaction region for the
514 // target.
515 //
516  G4int DspcT = (G4int) (AT*F) - Dabr;
517  if (DspcT <= 0) DspcT = 0;
518  else if (DspcT > AP-Dabr) DspcT = ((G4int) AT) - Dabr;
519 //
520 //
521 // Determine excitation energy associated with excess surface area of the
522 // target (EsT) and the excitation due to scattering of nucleons which are
523 // retained within the target (ExT). Add the total energy from the excited
524 // nucleus to the total energy of the secondaries.
525 //
526  if (fragmentT != NULL)
527  {
528  G4double EsT = theAbrasionGeometry->GetExcitationEnergyOfTarget();
529  G4double ExT = 0.0;
530  if (!excitationAbsorbedByProjectile)
531  ExT = GetNucleonInducedExcitation(rT, rP, r);
532  G4double xT = EsT + ExT;
533  if (xT > B*(AT-Dabr)) xT = B*(AT-Dabr);
534  G4LorentzVector lorentzVector = fragmentT->GetMomentum();
535  lorentzVector.setE(lorentzVector.e()+xT);
536  fragmentT->SetMomentum(lorentzVector);
537  TotalEPost += lorentzVector.e();
538  }
539 //
540 //
541 // Now determine the difference between the pre and post interaction
542 // energy - this will be used to determine the Lorentz boost if conservation
543 // of energy is to be imposed/attempted.
544 //
545  G4double deltaE = TotalEPre - TotalEPost;
546  if (deltaE > 0.0 && conserveEnergy)
547  {
548  G4double beta = std::sqrt(1.0 - EMassP*EMassP/G4Pow::GetInstance()->powN(deltaE+EMassP,2));
549  boost = boost / boost.mag() * beta;
550  }
551 //
552 //
553 // Now boost the secondaries from the projectile.
554 //
555  G4ThreeVector pBalance = pP.vect();
556  for (i=0; i<nSecP; i++)
557  {
559  GetParticle();
560  G4LorentzVector lorentzVector = dynamicP->Get4Momentum();
561  lorentzVector.boost(-boost);
562  dynamicP->Set4Momentum(lorentzVector);
563  pBalance -= lorentzVector.vect();
564  }
565 //
566 //
567 // Set the boost for the projectile prefragment. This is now based on the
568 // conservation of momentum. However, if the user selected momentum of the
569 // prefragment is not to be conserved this simply boosted to the velocity of the
570 // original projectile times the ratio of the unexcited to the excited mass
571 // of the prefragment (the excitation increases the effective mass of the
572 // prefragment, and therefore modifying the boost is an attempt to prevent
573 // the momentum of the prefragment being excessive).
574 //
575  if (fragmentP != NULL)
576  {
577  G4LorentzVector lorentzVector = fragmentP->GetMomentum();
578  G4double fragmentM = lorentzVector.m();
579  if (conserveMomentum)
580  fragmentP->SetMomentum
581  (G4LorentzVector(pBalance,std::sqrt(pBalance.mag2()+fragmentM*fragmentM+1.0*eV*eV)));
582  else
583  {
584  G4double fragmentGroundStateM = fragmentP->GetGroundStateMass();
585  fragmentP->SetMomentum(lorentzVector.boost(-boost * fragmentGroundStateM/fragmentM));
586  }
587  }
588 //
589 //
590 // Output information to user if verbose information requested.
591 //
592  if (verboseLevel >= 2)
593  {
594  G4cout <<G4endl;
595  G4cout <<"-----------------------------------" <<G4endl;
596  G4cout <<"Secondary nucleons from projectile:" <<G4endl;
597  G4cout <<"-----------------------------------" <<G4endl;
598  G4cout.precision(7);
599  for (i=0; i<nSecP; i++)
600  {
601  G4cout <<"Particle # " <<i <<G4endl;
604  G4cout <<"New nucleon (P) " <<dyn->GetDefinition()->GetParticleName()
605  <<" : " <<dyn->Get4Momentum()
606  <<G4endl;
607  }
608  G4cout <<"---------------------------" <<G4endl;
609  G4cout <<"The projectile prefragment:" <<G4endl;
610  G4cout <<"---------------------------" <<G4endl;
611  if (fragmentP != NULL)
612  G4cout <<*fragmentP <<G4endl;
613  else
614  G4cout <<"(No residual prefragment)" <<G4endl;
615  G4cout <<G4endl;
616  G4cout <<"-------------------------------" <<G4endl;
617  G4cout <<"Secondary nucleons from target:" <<G4endl;
618  G4cout <<"-------------------------------" <<G4endl;
619  G4cout.precision(7);
620  for (i=nSecP; i<nSec; i++)
621  {
622  G4cout <<"Particle # " <<i <<G4endl;
625  G4cout <<"New nucleon (T) " <<dyn->GetDefinition()->GetParticleName()
626  <<" : " <<dyn->Get4Momentum()
627  <<G4endl;
628  }
629  G4cout <<"-----------------------" <<G4endl;
630  G4cout <<"The target prefragment:" <<G4endl;
631  G4cout <<"-----------------------" <<G4endl;
632  if (fragmentT != NULL)
633  G4cout <<*fragmentT <<G4endl;
634  else
635  G4cout <<"(No residual prefragment)" <<G4endl;
636  }
637 //
638 //
639 // Now we can decay the nuclear fragments if present. The secondaries are
640 // collected and boosted as well. This is performed first for the projectile...
641 //
642  if (fragmentP !=NULL)
643  {
644  G4ReactionProductVector *products = NULL;
645  if (fragmentP->GetZ_asInt() != fragmentP->GetA_asInt())
646  products = theExcitationHandler->BreakItUp(*fragmentP);
647  else
648  products = theExcitationHandlerx->BreakItUp(*fragmentP);
649  delete fragmentP;
650  fragmentP = NULL;
651 
652  G4ReactionProductVector::iterator iter;
653  for (iter = products->begin(); iter != products->end(); ++iter)
654  {
655  G4DynamicParticle *secondary =
656  new G4DynamicParticle((*iter)->GetDefinition(),
657  (*iter)->GetTotalEnergy(), (*iter)->GetMomentum());
658  theParticleChange.AddSecondary (secondary); // Added MHM 20050118
659  G4String particleName = (*iter)->GetDefinition()->GetParticleName();
660  delete (*iter); // get rid of leftover particle def! // Added MHM 20050118
661  if (verboseLevel >= 2 && particleName.find("[",0) < particleName.size())
662  {
663  G4cout <<"------------------------" <<G4endl;
664  G4cout <<"The projectile fragment:" <<G4endl;
665  G4cout <<"------------------------" <<G4endl;
666  G4cout <<" fragmentP = " <<particleName
667  <<" Energy = " <<secondary->GetKineticEnergy()
668  <<G4endl;
669  }
670  }
671  delete products; // Added MHM 20050118
672  }
673 //
674 //
675 // Now decay the target nucleus - no boost is applied since in this
676 // approximation it is assumed that there is negligible momentum transfer from
677 // the projectile.
678 //
679  if (fragmentT != NULL)
680  {
681  G4ReactionProductVector *products = NULL;
682  if (fragmentT->GetZ_asInt() != fragmentT->GetA_asInt())
683  products = theExcitationHandler->BreakItUp(*fragmentT);
684  else
685  products = theExcitationHandlerx->BreakItUp(*fragmentT);
686  delete fragmentT;
687  fragmentT = NULL;
688 
689  G4ReactionProductVector::iterator iter;
690  for (iter = products->begin(); iter != products->end(); ++iter)
691  {
692  G4DynamicParticle *secondary =
693  new G4DynamicParticle((*iter)->GetDefinition(),
694  (*iter)->GetTotalEnergy(), (*iter)->GetMomentum());
695  theParticleChange.AddSecondary (secondary);
696  G4String particleName = (*iter)->GetDefinition()->GetParticleName();
697  delete (*iter); // get rid of leftover particle def! // Added MHM 20050118
698  if (verboseLevel >= 2 && particleName.find("[",0) < particleName.size())
699  {
700  G4cout <<"--------------------" <<G4endl;
701  G4cout <<"The target fragment:" <<G4endl;
702  G4cout <<"--------------------" <<G4endl;
703  G4cout <<" fragmentT = " <<particleName
704  <<" Energy = " <<secondary->GetKineticEnergy()
705  <<G4endl;
706  }
707  }
708  delete products; // Added MHM 20050118
709  }
710 
711  if (verboseLevel >= 2)
712  G4cout <<"########################################"
713  <<"########################################"
714  <<G4endl;
715 
716  delete theAbrasionGeometry;
717 
718  return &theParticleChange;
719 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double AtomicMass(const G4double A, const G4double Z) const
Definition: G4Nucleus.cc:254
G4long G4Poisson(G4double mean)
Definition: G4Poisson.hh:51
G4HadSecondary * GetSecondary(size_t i)
G4double GetKineticEnergy() const
long G4long
Definition: G4Types.hh:80
void DumpInfo(G4int mode=0) const
G4ParticleDefinition * GetDefinition() const
int G4int
Definition: G4Types.hh:78
G4ReactionProductVector * BreakItUp(const G4Fragment &theInitialState)
const G4String & GetParticleName() const
void SetStatusChange(G4HadFinalStateStatus aS)
std::vector< G4ReactionProduct * > G4ReactionProductVector
Hep3Vector vect() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
G4int GetA_asInt() const
Definition: G4Fragment.hh:266
const G4ParticleDefinition * GetDefinition() const
bool G4bool
Definition: G4Types.hh:79
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:307
HepLorentzVector & boost(double, double, double)
void SetMomentum(const G4LorentzVector &value)
Definition: G4Fragment.hh:312
G4double GetKineticEnergy() const
G4double GetEnergyDeposit()
Definition: G4Nucleus.hh:184
static constexpr double eV
Definition: G4SIunits.hh:215
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:288
const G4LorentzVector & Get4Momentum() const
G4LorentzVector Get4Momentum() const
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
void Set4Momentum(const G4LorentzVector &momentum)
void SetEnergyChange(G4double anEnergy)
Hep3Vector findBoostToCM() const
Hep3Vector unit() const
G4DynamicParticle * GetParticle()
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
double mag2() const
G4int GetZ_asInt() const
Definition: G4Fragment.hh:271
static constexpr double elm_coupling
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
static constexpr double fermi
Definition: G4SIunits.hh:103
G4double GetPDGCharge() const
double mag() const
void SetMomentumChange(const G4ThreeVector &aV)
G4int GetNumberOfSecondaries() const
G4double GetWilsonRadius(G4double A)
G4double GetTotalEnergy() const
CLHEP::HepLorentzVector G4LorentzVector

Here is the call graph for this function:

G4bool G4WilsonAbrasionModel::GetConserveMomentum ( )
inline

Definition at line 139 of file G4WilsonAbrasionModel.hh.

140  {return conserveMomentum;}
G4ExcitationHandler * G4WilsonAbrasionModel::GetExcitationHandler ( )
inline

Definition at line 124 of file G4WilsonAbrasionModel.hh.

125  {return theExcitationHandler;}
G4bool G4WilsonAbrasionModel::GetUseAblation ( )
inline

Definition at line 127 of file G4WilsonAbrasionModel.hh.

128  {return useAblation;}
void G4WilsonAbrasionModel::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 176 of file G4WilsonAbrasionModel.cc.

177 {
178  outFile << "G4WilsonAbrasionModel is a macroscopic treatment of\n"
179  << "nucleus-nucleus collisions using simple geometric arguments.\n"
180  << "The smaller projectile nucleus gouges out a part of the larger\n"
181  << "target nucleus, leaving a residual nucleus and a fireball\n"
182  << "region where the projectile and target intersect. The fireball"
183  << "is then treated as a highly excited nuclear fragment. This\n"
184  << "model is based on the NUCFRG2 model and is valid for all\n"
185  << "projectile energies between 70 MeV/n and 10.1 GeV/n. \n";
186 }
const G4WilsonAbrasionModel& G4WilsonAbrasionModel::operator= ( G4WilsonAbrasionModel right)
void G4WilsonAbrasionModel::SetConserveMomentum ( G4bool  conserveMomentum1)
inline

Definition at line 136 of file G4WilsonAbrasionModel.hh.

137  {conserveMomentum = conserveMomentum1;}
void G4WilsonAbrasionModel::SetExcitationHandler ( G4ExcitationHandler aExcitationHandler)
inline

Definition at line 121 of file G4WilsonAbrasionModel.hh.

122  {theExcitationHandler = aExcitationHandler;}
void G4WilsonAbrasionModel::SetUseAblation ( G4bool  useAblation1)

Definition at line 888 of file G4WilsonAbrasionModel.cc.

889 {
890  if (useAblation != useAblation1)
891  {
892  useAblation = useAblation1;
893  delete theExcitationHandler;
894  delete theExcitationHandlerx;
895  theExcitationHandler = new G4ExcitationHandler;
896  theExcitationHandlerx = new G4ExcitationHandler;
897  if (useAblation)
898  {
899  theAblation = new G4WilsonAblationModel;
900  theAblation->SetVerboseLevel(verboseLevel);
901  theExcitationHandler->SetEvaporation(theAblation);
902  theExcitationHandlerx->SetEvaporation(theAblation);
903  }
904  else
905  {
906  theAblation = NULL;
907  G4Evaporation * theEvaporation = new G4Evaporation;
908  G4FermiBreakUp * theFermiBreakUp = new G4FermiBreakUp;
909  G4StatMF * theMF = new G4StatMF;
910  theExcitationHandler->SetEvaporation(theEvaporation);
911  theExcitationHandler->SetFermiModel(theFermiBreakUp);
912  theExcitationHandler->SetMultiFragmentation(theMF);
913  theExcitationHandler->SetMaxAandZForFermiBreakUp(12, 6);
914  theExcitationHandler->SetMinEForMultiFrag(5.0*MeV);
915 
916  theEvaporation = new G4Evaporation;
917  theFermiBreakUp = new G4FermiBreakUp;
918  theExcitationHandlerx->SetEvaporation(theEvaporation);
919  theExcitationHandlerx->SetFermiModel(theFermiBreakUp);
920  theExcitationHandlerx->SetMaxAandZForFermiBreakUp(12, 6);
921  }
922  }
923  return;
924 }
void SetMinEForMultiFrag(G4double anE)
void SetFermiModel(G4VFermiBreakUp *ptr)
void SetMultiFragmentation(G4VMultiFragmentation *ptr)
void SetEvaporation(G4VEvaporation *ptr, G4bool isLocal=false)
void SetMaxAandZForFermiBreakUp(G4int anA, G4int aZ)
static constexpr double MeV
Definition: G4SIunits.hh:214

Here is the call graph for this function:

void G4WilsonAbrasionModel::SetVerboseLevel ( G4int  verboseLevel1)
inline

Definition at line 142 of file G4WilsonAbrasionModel.hh.

143 {
144  verboseLevel = verboseLevel1;
145  if (useAblation) theAblation->SetVerboseLevel(verboseLevel);
146 }

Here is the call graph for this function:


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