Geant4  10.00.p02
G4GeneratorPrecompoundInterface.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4GeneratorPrecompoundInterface.cc 66812 2013-01-12 16:06:46Z gcosmo $
27 //
28 // -----------------------------------------------------------------------------
29 // GEANT 4 class file
30 //
31 // History: first implementation
32 // HPW, 10DEC 98, the decay part originally written by Gunter Folger
33 // in his FTF-test-program.
34 //
35 // M.Kelsey, 28 Jul 2011 -- Replace loop to decay input secondaries
36 // with new utility class, simplify cleanup loops
37 // -----------------------------------------------------------------------------
38 
39 #include <algorithm>
40 #include <vector>
41 
43 #include "G4PhysicalConstants.hh"
44 #include "G4SystemOfUnits.hh"
46 #include "G4KineticTrackVector.hh"
47 #include "G4Proton.hh"
48 #include "G4Neutron.hh"
49 
50 #include "G4Deuteron.hh"
51 #include "G4Triton.hh"
52 #include "G4He3.hh"
53 #include "G4Alpha.hh"
54 
55 #include "G4V3DNucleus.hh"
56 #include "G4Nucleon.hh"
57 
58 #include "G4AntiProton.hh"
59 #include "G4AntiNeutron.hh"
60 #include "G4AntiDeuteron.hh"
61 #include "G4AntiTriton.hh"
62 #include "G4AntiHe3.hh"
63 #include "G4AntiAlpha.hh"
64 
65 #include "G4FragmentVector.hh"
66 #include "G4ReactionProduct.hh"
68 #include "G4PreCompoundModel.hh"
69 #include "G4ExcitationHandler.hh"
70 #include "G4DecayKineticTracks.hh"
72 
74 : CaptureThreshold(10*MeV)
75 {
78 
81  He3 =G4He3::He3();
83 
86 
91 
92  if(preModel) { SetDeExcitation(preModel); }
93  else {
94  G4HadronicInteraction* hadi =
96  G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(hadi);
97  if(!pre) { pre = new G4PreCompoundModel(); }
98  SetDeExcitation(pre);
99  }
100 }
101 
103 {
104 }
105 
106 //---------------------------------------------------------------------
107 // choose to calculate excitation energy from energy balance
108 #define exactExcitationEnergy
109 //#define debugPrecoInt
110 //#define G4GPI_debug_excitation
111 
113 Propagate(G4KineticTrackVector* theSecondaries, G4V3DNucleus* theNucleus)
114 {
115  G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
116 
117  // decay the strong resonances
118  G4DecayKineticTracks decay(theSecondaries);
119 
120  // prepare the fragment
121  G4int anA=theNucleus->GetMassNumber();
122  G4int aZ=theNucleus->GetCharge();
123  G4int numberOfEx = 0;
124  G4int numberOfCh = 0;
125  G4int numberOfHoles = 0;
126  G4double exEnergy = 0.0;
127  G4double R = theNucleus->GetNuclearRadius();
128  G4ThreeVector exciton3Momentum(0.,0.,0.);
129  G4ThreeVector captured3Momentum(0.,0.,0.);
130  G4ThreeVector wounded3Momentum(0.,0.,0.);
131 
132  // loop over secondaries
133 #ifdef exactExcitationEnergy
134  G4LorentzVector secondary4Momemtum(0,0,0,0);
135 #endif
136  G4KineticTrackVector::iterator iter;
137  for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
138  {
139  G4ParticleDefinition* part = (*iter)->GetDefinition();
140  G4double e = (*iter)->Get4Momentum().e();
141  G4double mass = (*iter)->Get4Momentum().mag();
142  G4ThreeVector mom = (*iter)->Get4Momentum().vect();
143  if((part != proton && part != neutron) ||
144  (e > mass + CaptureThreshold) ||
145  ((*iter)->GetPosition().mag() > R)) {
146  G4ReactionProduct * theNew = new G4ReactionProduct(part);
147  theNew->SetMomentum(mom);
148  theNew->SetTotalEnergy(e);
149  theTotalResult->push_back(theNew);
150 #ifdef exactExcitationEnergy
151  secondary4Momemtum += (*iter)->Get4Momentum();
152 #endif
153  } else {
154  // within the nucleus, neutron or proton
155  // now calculate A, Z of the fragment, momentum, number of exciton states
156  ++anA;
157  ++numberOfEx;
158  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
159  aZ += Z;
160  numberOfCh += Z;
161  captured3Momentum += mom;
162  exEnergy += (e - mass);
163  }
164  delete (*iter);
165  }
166  delete theSecondaries;
167 
168  // loop over wounded nucleus
169  G4Nucleon * theCurrentNucleon =
170  theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
171  while(theCurrentNucleon) {
172  if(theCurrentNucleon->AreYouHit()) {
173  ++numberOfHoles;
174  ++numberOfEx;
175  --anA;
176  aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/eplus + 0.1);
177  wounded3Momentum += theCurrentNucleon->Get4Momentum().vect();
178  //G4cout << "hit nucleon " << theCurrentNucleon->Get4Momentum() << G4endl;
179  exEnergy += theCurrentNucleon->GetBindingEnergy();
180  }
181  theCurrentNucleon = theNucleus->GetNextNucleon();
182  }
183  exciton3Momentum = captured3Momentum - wounded3Momentum;
184 
185  if(anA == 0) return theTotalResult;
186 
187  if(anA >= aZ)
188  {
190 
191 #ifdef exactExcitationEnergy
192  // recalculate exEnergy from Energy balance....
193  const G4HadProjectile * primary = GetPrimaryProjectile();
194  G4double Einitial= primary->Get4Momentum().e()
196  theNucleus->GetCharge());
197 // Uzhi G4double Efinal = fMass + secondary4Momemtum.e();
198  G4double Efinal = std::sqrt(exciton3Momentum.mag2() + fMass*fMass)
199  + secondary4Momemtum.e();
200  if ( (Einitial - Efinal) > 0 ) {
201  // G4cout << "G4GPI::Propagate() : positive exact excitation Energy "
202  // << (Einitial - Efinal)/MeV << " MeV, exciton estimate "
203  // << exEnergy/MeV << " MeV" << G4endl;
204 
205 // exEnergy=Einitial - Efinal;
206  G4LorentzVector PrimMom=primary->Get4Momentum(); PrimMom.setE(Einitial);
207 
208  exEnergy=(PrimMom - secondary4Momemtum).mag() - fMass;
209  }
210  else {
211  // G4cout << "G4GeneratorPrecompoundInterface::Propagate() : "
212  // << "negative exact excitation Energy "
213  // << (Einitial - Efinal)/MeV
214  // << " MeV, setting excitation to 0 MeV" << G4endl;
215  exEnergy=0.;
216  }
217 #endif
218 
219  if(exEnergy < 0.) exEnergy=0.; // Uzhi 11 Dec. 2012
220 
221  fMass += exEnergy;
222 
223  G4ThreeVector balance=primary->Get4Momentum().vect() -
224  secondary4Momemtum.vect() - exciton3Momentum;
225 
226 #ifdef G4GPI_debug_excitation
227  G4cout << "momentum balance" << balance
228  << " value " << balance.mag() <<G4endl
229  << "primary "<< primary->Get4Momentum() <<G4endl
230  << "secondary "<< secondary4Momemtum <<G4endl
231  << "captured "<< captured3Momentum <<G4endl
232  << "wounded "<< wounded3Momentum <<G4endl
233  << "exciton "<< exciton3Momentum <<G4endl
234  << "second + exciton"
235  << secondary4Momemtum.vect() + exciton3Momentum << G4endl;
236 #endif
237 //#ifdef exactExcitationEnergy
238 // G4LorentzVector exciton4Momentum(exciton3Momentum, fMass);
239 // G4LorentzVector exciton4Momentum(exciton3Momentum,
240 // std::sqrt(exciton3Momentum.mag2() + fMass*fMass));
241 //#else
242  G4LorentzVector exciton4Momentum(exciton3Momentum,
243  std::sqrt(exciton3Momentum.mag2() + fMass*fMass));
244 //#endif
245 //G4cout<<"exciton4Momentum "<<exciton4Momentum<<G4endl;
246  // Need to de-excite the remnant nucleus only if excitation energy > 0.
247  G4Fragment anInitialState(anA, aZ, exciton4Momentum);
248  anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
249  anInitialState.SetNumberOfCharged(numberOfCh);
250  anInitialState.SetNumberOfHoles(numberOfHoles);
251 
252  G4ReactionProductVector * aPrecoResult =
253  theDeExcitation->DeExcite(anInitialState);
254  // fill pre-compound part into the result, and return
255  theTotalResult->insert(theTotalResult->end(),aPrecoResult->begin(),
256  aPrecoResult->end() );
257  delete aPrecoResult;
258  }
259 
260  return theTotalResult;
261 }
262 
265 {
266  G4cout << "G4GeneratorPrecompoundInterface: ApplyYourself interface called stand-allone."
267  << G4endl;
268  G4cout << "This class is only a mediator between generator and precompound"<<G4endl;
269  G4cout << "Please remove from your physics list."<<G4endl;
270  throw G4HadronicException(__FILE__, __LINE__, "SEVERE: G4GeneratorPrecompoundInterface model interface called stand-allone.");
271  return new G4HadFinalState;
272 }
274 {
275  outFile << "G4GeneratorPrecompoundInterface interfaces a high\n"
276  << "energy model through the wounded nucleus to precompound de-excition.\n"
277  << "Low energy protons and neutron present among secondaries produced by \n"
278  << "the high energy generator and within the nucleus are captured. The wounded\n"
279  << "nucleus and the captured particles form an excited nuclear fragment. This\n"
280  << "fragment is passed to the Geant4 pre-compound model for de-excitation.\n"
281  << "Nuclear de-excitation:\n";
282  // preco
283 
284 }
285 
286 // Uzhi Nov. 2012 ------------------------------------------------
289  G4V3DNucleus* theProjectileNucleus)
290 {
291 #ifdef debugPrecoInt
292  G4cout<<"G4GeneratorPrecompoundInterface::PropagateNuclNucl "<<G4endl;
293 #endif
294 
295  G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
296 
297  // prepare the target residual
298  G4int anA=theNucleus->GetMassNumber();
299  G4int aZ=theNucleus->GetCharge();
300  G4int numberOfEx = 0;
301  G4int numberOfCh = 0;
302  G4int numberOfHoles = 0;
303  G4double exEnergy = 0.0;
304  G4double R = theNucleus->GetNuclearRadius();
305  G4LorentzVector Target4Momentum(0,0,0,0);
306 
307 #ifdef debugPrecoInt
308  G4cout<<"Target A Z "<<anA<<" "<<aZ<<G4endl;
309 #endif
310 
311  // loop over wounded target nucleus
312  G4Nucleon * theCurrentNucleon =
313  theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
314  while(theCurrentNucleon) {
315  if(theCurrentNucleon->AreYouHit()) {
316  ++numberOfHoles;
317  ++numberOfEx;
318  --anA;
319  aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
320  eplus + 0.1);
321  exEnergy += theCurrentNucleon->GetBindingEnergy();
322  Target4Momentum -=theCurrentNucleon->Get4Momentum();
323  }
324  theCurrentNucleon = theNucleus->GetNextNucleon();
325  }
326 
327 #ifdef debugPrecoInt
328  G4cout<<"Residual Target A Z E* 4mom "<<anA<<" "<<aZ<<" "<<exEnergy<<" "
329  <<Target4Momentum<<G4endl;
330 #endif
331 
332  // prepare the projectile residual
333 #ifdef debugPrecoInt
334  G4cout<<"Primary BaryonNumber "
336 #endif
337 
338  G4bool ProjectileIsAntiNucleus=
340 
341  G4ThreeVector bst = GetPrimaryProjectile()->Get4Momentum().boostVector();
342 
343  G4int anAb=theProjectileNucleus->GetMassNumber();
344  G4int aZb=theProjectileNucleus->GetCharge();
345  G4int numberOfExB = 0;
346  G4int numberOfChB = 0;
347  G4int numberOfHolesB = 0;
348  G4double exEnergyB = 0.0;
349  G4double Rb = theProjectileNucleus->GetNuclearRadius();
350  G4LorentzVector Projectile4Momentum(0,0,0,0);
351 
352 #ifdef debugPrecoInt
353  G4cout<<"Projectile A Z "<<anAb<<" "<<aZb<<G4endl;
354 #endif
355 
356  // loop over wounded projectile nucleus
357  theCurrentNucleon =
358  theProjectileNucleus->StartLoop() ? theProjectileNucleus->GetNextNucleon() : 0;
359  while(theCurrentNucleon) {
360  if(theCurrentNucleon->AreYouHit()) {
361  ++numberOfHolesB;
362  ++numberOfExB;
363  --anAb;
364  if(!ProjectileIsAntiNucleus)
365  {
366  aZb -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
367  eplus + 0.1);
368  } else
369  {
370  aZb += G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
371  eplus - 0.1);
372  }
373  exEnergyB += theCurrentNucleon->GetBindingEnergy();
374  Projectile4Momentum -=theCurrentNucleon->Get4Momentum();
375  }
376  theCurrentNucleon = theProjectileNucleus->GetNextNucleon();
377  }
378 
379  G4bool ExistTargetRemnant = G4double (numberOfHoles) <
380  0.3* G4double (numberOfHoles + anA);
381  G4bool ExistProjectileRemnant= G4double (numberOfHolesB) <
382  0.3*G4double (numberOfHolesB + anAb);
383 
384 #ifdef debugPrecoInt
385  G4cout<<"Projectile residual A Z E* 4mom "<<anAb<<" "<<aZb<<" "<<exEnergyB<<" "
386  <<Projectile4Momentum<<G4endl;
387  G4cout<<" ExistTargetRemnant ExistProjectileRemnant "
388  <<ExistTargetRemnant<<" "<< ExistProjectileRemnant<<G4endl;
389 #endif
390 //-----------------------------------------------------------------------------
391  // decay the strong resonances
392  G4DecayKineticTracks decay(theSecondaries);
393 
394 #ifdef debugPrecoInt
395  G4LorentzVector secondary4Momemtum(0,0,0,0);
396  G4int SecondrNum(0);
397 #endif
398 
399  // loop over secondaries
400  G4KineticTrackVector::iterator iter;
401  for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
402  {
403  G4ParticleDefinition* part = (*iter)->GetDefinition();
404  G4LorentzVector aTrack4Momentum=(*iter)->Get4Momentum();
405 
406  if( part != proton && part != neutron &&
407  (part != ANTIproton && ProjectileIsAntiNucleus) &&
408  (part != ANTIneutron && ProjectileIsAntiNucleus) )
409  {
410  G4ReactionProduct * theNew = new G4ReactionProduct(part);
411  theNew->SetMomentum(aTrack4Momentum.vect());
412  theNew->SetTotalEnergy(aTrack4Momentum.e());
413  theTotalResult->push_back(theNew);
414 #ifdef debugPrecoInt
415  SecondrNum++;
416  secondary4Momemtum += (*iter)->Get4Momentum();
417  G4cout<<"Secondary "<<SecondrNum<<" "
418  <<theNew->GetDefinition()->GetParticleName()<<" "
419  <<secondary4Momemtum<<G4endl;
420 #endif
421  delete (*iter);
422  continue;
423  }
424 
425  G4bool CanBeCapturedByTarget = false;
426  if( part == proton || part == neutron)
427  {
428  CanBeCapturedByTarget = ExistTargetRemnant &&
430  (aTrack4Momentum + Target4Momentum).mag() -
431  aTrack4Momentum.mag() - Target4Momentum.mag()) &&
432  ((*iter)->GetPosition().mag() < R);
433  }
434 // ---------------------------
435  G4LorentzVector Position((*iter)->GetPosition(),
436  (*iter)->GetFormationTime());
437  Position.boost(bst);
438 
439  G4bool CanBeCapturedByProjectile = false;
440 
441  if( !ProjectileIsAntiNucleus &&
442  ( part == proton || part == neutron))
443  {
444  CanBeCapturedByProjectile = ExistProjectileRemnant &&
446  (aTrack4Momentum + Projectile4Momentum).mag() -
447  aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
448  (Position.vect().mag() < Rb);
449  }
450 
451  if( ProjectileIsAntiNucleus &&
452  ( part == ANTIproton || part == ANTIneutron))
453  {
454  CanBeCapturedByProjectile = ExistProjectileRemnant &&
456  (aTrack4Momentum + Projectile4Momentum).mag() -
457  aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
458  (Position.vect().mag() < Rb);
459  }
460 
461  if(CanBeCapturedByTarget && CanBeCapturedByProjectile)
462  {
463  if(G4UniformRand() < 0.5)
464  { CanBeCapturedByTarget = true; CanBeCapturedByProjectile = false;}
465  else
466  { CanBeCapturedByTarget = false; CanBeCapturedByProjectile = true;}
467  }
468 
469  if(CanBeCapturedByTarget)
470  {
471  // within the target nucleus, neutron or proton
472  // now calculate A, Z of the fragment, momentum,
473  // number of exciton states
474 #ifdef debugPrecoInt
475  G4cout<<"Track is CapturedByTarget "<<" "
476  <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
477 #endif
478  ++anA;
479  ++numberOfEx;
480  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
481  aZ += Z;
482  numberOfCh += Z;
483  Target4Momentum +=aTrack4Momentum;
484  delete (*iter);
485  } else if(CanBeCapturedByProjectile)
486  {
487  // within the projectile nucleus, neutron or proton
488  // now calculate A, Z of the fragment, momentum,
489  // number of exciton states
490 #ifdef debugPrecoInt
491  G4cout<<"Track is CapturedByProjectile"<<" "
492  <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
493 #endif
494  ++anAb;
495  ++numberOfExB;
496  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
497  if( ProjectileIsAntiNucleus ) Z=-Z;
498  aZb += Z;
499  numberOfChB += Z;
500  Projectile4Momentum +=aTrack4Momentum;
501  delete (*iter);
502  } else
503  { // the track is not captured
504  G4ReactionProduct * theNew = new G4ReactionProduct(part);
505  theNew->SetMomentum(aTrack4Momentum.vect());
506  theNew->SetTotalEnergy(aTrack4Momentum.e());
507  theTotalResult->push_back(theNew);
508 
509 #ifdef debugPrecoInt
510  SecondrNum++;
511  secondary4Momemtum += (*iter)->Get4Momentum();
512  G4cout<<"Secondary "<<SecondrNum<<" "
513  <<theNew->GetDefinition()->GetParticleName()<<" "
514  <<secondary4Momemtum<<G4endl;
515 #endif
516  delete (*iter);
517  continue;
518  }
519  }
520  delete theSecondaries;
521 //-----------------------------------------------------
522 
523 #ifdef debugPrecoInt
524  G4cout<<"Final target residual A Z E* 4mom "<<anA<<" "<<aZ<<" "
525  <<exEnergy<<" "<<Target4Momentum<<G4endl;
526 #endif
527 
528  if(0!=anA )
529  {
531 
532  if((anA == theNucleus->GetMassNumber()) && (exEnergy <= 0.))
533  {Target4Momentum.setE(fMass);}
534 
535  G4double RemnMass=Target4Momentum.mag();
536  if(RemnMass < fMass)
537  {
538  RemnMass=fMass + exEnergy;
539  Target4Momentum.setE(std::sqrt(Target4Momentum.vect().mag2() +
540  RemnMass*RemnMass));
541  } else
542  { exEnergy=RemnMass-fMass;}
543 
544  if( exEnergy < 0.) exEnergy=0.;
545 
546  // Need to de-excite the remnant nucleus
547  G4Fragment anInitialState(anA, aZ, Target4Momentum);
548  anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
549  anInitialState.SetNumberOfCharged(numberOfCh);
550  anInitialState.SetNumberOfHoles(numberOfHoles);
551 
552  G4ReactionProductVector * aPrecoResult =
553  theDeExcitation->DeExcite(anInitialState);
554 
555  // fill pre-compound part into the result, and return
556  for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
557  {
558  theTotalResult->push_back(aPrecoResult->operator[](ll));
559 #ifdef debugPrecoInt
560  G4cout<<"Tr frag "<<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()
561  <<" "<<aPrecoResult->operator[](ll)->GetMomentum()<<G4endl;
562 #endif
563  }
564  delete aPrecoResult;
565  }
566 
567 //-----------------------------------------------------
568  if((anAb == theProjectileNucleus->GetMassNumber())&& (exEnergyB <= 0.))
569  {Projectile4Momentum = GetPrimaryProjectile()->Get4Momentum();}
570 
571 #ifdef debugPrecoInt
572  G4cout<<"Final projectile residual A Z E* Pmom "<<anAb<<" "<<aZb<<" "
573  <<exEnergyB<<" "<<Projectile4Momentum<<G4endl;
574 #endif
575 
576  if(0!=anAb)
577  {
578  G4double fMass = G4NucleiProperties::GetNuclearMass(anAb, aZb);
579  G4double RemnMass=Projectile4Momentum.mag();
580 
581  if(RemnMass < fMass)
582  {
583  RemnMass=fMass + exEnergyB;
584  Projectile4Momentum.setE(std::sqrt(Projectile4Momentum.vect().mag2() +
585  RemnMass*RemnMass));
586  } else
587  { exEnergyB=RemnMass-fMass;}
588 
589  if( exEnergyB < 0.) exEnergyB=0.;
590 
591  // Need to de-excite the remnant nucleus
592  G4Fragment anInitialState(anAb, aZb, Projectile4Momentum);
593  anInitialState.SetNumberOfParticles(numberOfExB-numberOfHolesB);
594  anInitialState.SetNumberOfCharged(numberOfChB);
595  anInitialState.SetNumberOfHoles(numberOfHolesB);
596 
597  G4ReactionProductVector * aPrecoResult =
598  theDeExcitation->DeExcite(anInitialState);
599 
600  // fill pre-compound part into the result, and return
601  for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
602  {
603  if(ProjectileIsAntiNucleus)
604  {
605 
606 #ifdef debugPrecoInt
607  G4cout<<"aPrecoRes "<<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()
608  <<" "<<aPrecoResult->operator[](ll)->GetMomentum()
609  <<" "<<aPrecoResult->operator[](ll)->GetTotalEnergy()
610  <<" "<<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
611 #endif
612 
613  G4ParticleDefinition * aFragment=aPrecoResult->operator[](ll)->GetDefinition();
614  G4ParticleDefinition * LastFragment=aFragment;
615  if (aFragment == proton) {LastFragment=G4AntiProton::AntiProtonDefinition();}
616  else if(aFragment == neutron) {LastFragment=G4AntiNeutron::AntiNeutronDefinition();}
617  else if(aFragment == deuteron){LastFragment=G4AntiDeuteron::AntiDeuteronDefinition();}
618  else if(aFragment == triton) {LastFragment=G4AntiTriton::AntiTritonDefinition();}
619  else if(aFragment == He3) {LastFragment=G4AntiHe3::AntiHe3Definition();}
620  else if(aFragment == He4) {LastFragment=G4AntiAlpha::AntiAlphaDefinition();}
621  else {}
622 
623  aPrecoResult->operator[](ll)->SetDefinitionAndUpdateE(LastFragment);
624  }
625 
626 #ifdef debugPrecoInt
627  G4cout<<"aPrecoResA "<<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()
628  <<" "<<aPrecoResult->operator[](ll)->GetMomentum()
629  <<" "<<aPrecoResult->operator[](ll)->GetTotalEnergy()
630  <<" "<<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
631 #endif
632  theTotalResult->push_back(aPrecoResult->operator[](ll));
633  }
634 
635  delete aPrecoResult;
636  }
637 
638  return theTotalResult;
639 }
640 
641 // Uzhi Nov. 2012 ------------------------------------------------
642 
static G4AntiTriton * AntiTritonDefinition()
Definition: G4AntiTriton.cc:89
static G4AntiHe3 * AntiHe3()
Definition: G4AntiHe3.cc:94
static const double MeV
Definition: G4SIunits.hh:193
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4GeneratorPrecompoundInterface(G4VPreCompoundModel *p=0)
virtual G4int GetCharge()=0
CLHEP::Hep3Vector G4ThreeVector
std::ofstream outFile
Definition: GammaRayTel.cc:68
const G4HadProjectile * GetPrimaryProjectile() const
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)=0
virtual G4double GetNuclearRadius()=0
virtual G4ReactionProductVector * Propagate(G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
virtual G4bool StartLoop()=0
static G4AntiDeuteron * AntiDeuteron()
void SetMomentum(const G4double x, const G4double y, const G4double z)
virtual G4int GetMassNumber()=0
virtual void PropagateModelDescription(std::ostream &) const
void SetNumberOfHoles(G4int valueTot, G4int valueP=0)
Definition: G4Fragment.hh:355
static G4AntiDeuteron * AntiDeuteronDefinition()
virtual const G4LorentzVector & Get4Momentum() const
Definition: G4Nucleon.hh:72
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
static G4AntiAlpha * AntiAlpha()
Definition: G4AntiAlpha.cc:89
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:88
int G4int
Definition: G4Types.hh:78
static G4AntiNeutron * AntiNeutronDefinition()
const G4String & GetParticleName() const
G4ParticleDefinition * GetDefinition() const
std::vector< G4ReactionProduct * > G4ReactionProductVector
#define G4UniformRand()
Definition: Randomize.hh:87
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
const G4ParticleDefinition * GetDefinition() const
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
bool G4bool
Definition: G4Types.hh:79
G4bool AreYouHit() const
Definition: G4Nucleon.hh:97
void SetTotalEnergy(const G4double en)
static G4Triton * Triton()
Definition: G4Triton.cc:95
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
void SetNumberOfParticles(G4int value)
Definition: G4Fragment.hh:364
virtual G4ReactionProductVector * PropagateNuclNucl(G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4V3DNucleus *theProjectileNucleus)
const G4LorentzVector & Get4Momentum() const
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4AntiHe3 * AntiHe3Definition()
Definition: G4AntiHe3.cc:89
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
void SetDeExcitation(G4VPreCompoundModel *ptr)
static G4AntiAlpha * AntiAlphaDefinition()
Definition: G4AntiAlpha.cc:84
#define G4endl
Definition: G4ios.hh:61
virtual G4ParticleDefinition * GetDefinition() const
Definition: G4Nucleon.hh:85
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
virtual G4Nucleon * GetNextNucleon()=0
void SetNumberOfCharged(G4int value)
Definition: G4Fragment.hh:369
double G4double
Definition: G4Types.hh:76
static G4AntiTriton * AntiTriton()
Definition: G4AntiTriton.cc:94
static const double eplus
Definition: G4SIunits.hh:178
G4double GetPDGCharge() const
static G4He3 * He3()
Definition: G4He3.cc:94
G4double GetBindingEnergy() const
Definition: G4Nucleon.hh:75
static G4AntiNeutron * AntiNeutron()
CLHEP::HepLorentzVector G4LorentzVector