Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4NeutronRadCapture.cc
이 파일의 문서화 페이지로 가기
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: G4NeutronRadCapture.cc 88191 2015-02-02 17:27:37Z gcosmo $
27 //
28 //
29 // Physics model class G4NeutronRadCapture
30 // Created: 31 August 2009
31 // Author V.Ivanchenko
32 //
33 // Modified:
34 // 09.09.2010 V.Ivanchenko added usage of G4PhotonEvaporation
35 //
36 
37 #include "G4NeutronRadCapture.hh"
38 #include "G4SystemOfUnits.hh"
39 #include "G4ParticleDefinition.hh"
40 #include "G4Fragment.hh"
41 #include "G4FragmentVector.hh"
42 #include "G4NucleiProperties.hh"
43 #include "G4VEvaporationChannel.hh"
44 #include "G4PhotonEvaporation.hh"
45 #include "G4DynamicParticle.hh"
46 #include "G4ParticleTable.hh"
47 #include "G4IonTable.hh"
48 #include "G4Electron.hh"
49 #include "G4Deuteron.hh"
50 #include "G4Triton.hh"
51 #include "G4He3.hh"
52 #include "G4Alpha.hh"
53 #include "G4RandomDirection.hh"
54 
56  : G4HadronicInteraction("nRadCapture"),
57  photonEvaporation(nullptr),lab4mom(0.,0.,0.,0.)
58 {
61  SetMinEnergy( 0.0*CLHEP::GeV );
62  SetMaxEnergy( 100.*CLHEP::TeV );
63 
65  icID = -1;
66 
68 }
69 
71 {
72  delete photonEvaporation;
73 }
74 
76 {
77  if(photonEvaporation != nullptr) { return; }
78  G4DeexPrecoParameters* param =
81  icID = param->GetInternalConversionID();
82 
86 }
87 
89  const G4HadProjectile& aTrack, G4Nucleus& theNucleus)
90 {
93 
94  G4int A = theNucleus.GetA_asInt();
95  G4int Z = theNucleus.GetZ_asInt();
96 
97  G4double time = aTrack.GetGlobalTime();
98 
99  // Create initial state
101  lab4mom += aTrack.Get4Momentum();
102 
103  G4double M = lab4mom.mag();
104  ++A;
106  //G4cout << "Capture start: Z= " << Z << " A= " << A
107  // << " LabM= " << M << " Mcompound= " << mass << G4endl;
108 
109  // simplified method of 1 gamma emission
110  if(A <= 4) {
111 
113 
114  if(M - mass <= lowestEnergyLimit) {
115  return &theParticleChange;
116  }
117 
118  if (verboseLevel > 1) {
119  G4cout << "G4NeutronRadCapture::DoIt: Eini(MeV)="
120  << aTrack.GetKineticEnergy()/MeV << " Eexc(MeV)= "
121  << (M - mass)/MeV
122  << " Z= " << Z << " A= " << A << G4endl;
123  }
124  G4double e1 = (M - mass)*(M + mass)/(2*M);
125  G4LorentzVector lv2(e1*G4RandomDirection(),e1);
126  lv2.boost(bst);
127  G4HadSecondary* news =
129  news->SetTime(time);
131  delete news;
132 
133  const G4ParticleDefinition* theDef = 0;
134 
135  lab4mom -= lv2;
136  if (Z == 1 && A == 2) {theDef = G4Deuteron::Deuteron();}
137  else if (Z == 1 && A == 3) {theDef = G4Triton::Triton();}
138  else if (Z == 2 && A == 3) {theDef = G4He3::He3();}
139  else if (Z == 2 && A == 4) {theDef = G4Alpha::Alpha();}
140  else { theDef = theTableOfIons->GetIon(Z,A,0.0,noFloat,0); }
141 
142  if (verboseLevel > 1) {
143  G4cout << "Gamma 4-mom: " << lv2 << " "
144  << theDef->GetParticleName() << " " << lab4mom << G4endl;
145  }
146  if(theDef) {
147  news = new G4HadSecondary(new G4DynamicParticle(theDef, lab4mom));
148  news->SetTime(time);
150  delete news;
151  }
152 
153  // Use photon evaporation
154  } else {
155 
156  // protection against wrong kinematic
157  if(M < mass) {
158  G4double etot = std::max(mass, lab4mom.e());
159  G4double ptot = std::sqrt((etot - mass)*(etot + mass));
160  G4ThreeVector v = lab4mom.vect().unit();
161  lab4mom.set(v.x()*ptot,v.y()*ptot,v.z()*ptot,etot);
162  }
163 
164  G4Fragment* aFragment = new G4Fragment(A, Z, lab4mom);
165 
166  if (verboseLevel > 1) {
167  G4cout << "G4NeutronRadCapture::ApplyYourself initial G4Fragmet:"
168  << G4endl;
169  G4cout << aFragment << G4endl;
170  }
171 
172  //
173  // Sample final state
174  //
176  if(!fv) { fv = new G4FragmentVector(); }
177  fv->push_back(aFragment);
178  size_t n = fv->size();
179 
180  if (verboseLevel > 1) {
181  G4cout << "G4NeutronRadCapture: " << n << " final particle icID= " << icID << G4endl;
182  }
183  for(size_t i=0; i<n; ++i) {
184 
185  G4Fragment* f = (*fv)[i];
186  G4double etot = f->GetMomentum().e();
187 
188  Z = f->GetZ_asInt();
189  A = f->GetA_asInt();
190 
191  const G4ParticleDefinition* theDef;
192  if(0 == Z && 0 == A) {theDef = f->GetParticleDefinition();}
193  else if (Z == 1 && A == 2) {theDef = G4Deuteron::Deuteron();}
194  else if (Z == 1 && A == 3) {theDef = G4Triton::Triton();}
195  else if (Z == 2 && A == 3) {theDef = G4He3::He3();}
196  else if (Z == 2 && A == 4) {theDef = G4Alpha::Alpha();}
197  else {
198  G4double eexc = f->GetExcitationEnergy();
199  if(eexc <= minExcitation) { eexc = 0.0; }
200  theDef = theTableOfIons->GetIon(Z, A, eexc, noFloat, 0);
201  /*
202  G4cout << "### NC Find ion Z= " << Z << " A= " << A
203  << " Eexc(MeV)= " << eexc/MeV << " "
204  << theDef << G4endl;
205  */
206  }
207  G4double ekin = std::max(0.0,etot - theDef->GetPDGMass());
208  if (verboseLevel > 1) {
209  G4cout << i << ". " << theDef->GetParticleName()
210  << " Ekin(MeV)= " << etot/MeV
211  << " p: " << f->GetMomentum().vect()
212  << G4endl;
213  }
214  G4HadSecondary* news = new G4HadSecondary(
215  new G4DynamicParticle(theDef,
216  f->GetMomentum().vect().unit(),
217  ekin));
218  G4double timeF = f->GetCreationTime();
219  if(timeF < 0.0) { timeF = 0.0; }
220  news->SetTime(time + timeF);
221  if(theDef == electron) { news->SetCreatorModelType(icID); }
223  delete news;
224  delete f;
225  }
226  delete fv;
227  }
228  //G4cout << "Capture done" << G4endl;
229  return &theParticleChange;
230 }
231 
static G4He3 * He3()
Definition: G4He3.cc:94
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double GetCreationTime() const
Definition: G4Fragment.hh:441
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:300
void SetMinEnergy(G4double anEnergy)
void SetMaxEnergy(const G4double anEnergy)
static G4ParticleTable * GetParticleTable()
static constexpr double MeV
Definition: G4SIunits.hh:214
G4DeexPrecoParameters * GetParameters()
#define G4endl
Definition: G4ios.hh:61
static G4NuclearLevelData * GetInstance()
G4int GetInternalConversionID() const
const G4String & GetParticleName() const
G4int GetA_asInt() const
Definition: G4Fragment.hh:259
G4IonTable * GetIonTable() const
double z() const
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
G4ThreeVector G4RandomDirection()
static constexpr double TeV
G4double GetMinExcitation() const
G4double GetPDGMass() const
G4FragmentVector * BreakUpFragment(G4Fragment *theNucleus)
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
Float_t Z
double G4double
Definition: G4Types.hh:76
G4double GetGlobalTime() const
virtual void InitialiseModel() final
void SetCreatorModelType(G4int idx)
double A(double temperature)
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4double GetNuclearMass(const G4double A, const G4double Z)
static constexpr double eV
void set(double x, double y, double z, double t)
static G4Electron * Electron()
Definition: G4Electron.cc:94
#define noFloat
Definition: G4Ions.hh:118
G4int GetZ_asInt() const
Definition: G4Fragment.hh:264
const G4LorentzVector & Get4Momentum() const
Hep3Vector unit() const
virtual void SetICM(G4bool)
G4double GetKineticEnergy() const
double mag() const
int G4int
Definition: G4Types.hh:78
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
static G4Triton * Triton()
Definition: G4Triton.cc:95
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus) final
G4VEvaporationChannel * photonEvaporation
G4GLOB_DLL std::ostream G4cout
double x() const
Hep3Vector boostVector() const
const G4ParticleDefinition * GetParticleDefinition() const
Definition: G4Fragment.hh:431
Hep3Vector vect() const
Char_t n[5]
void SetTime(G4double aT)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:504
double y() const
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:276
static constexpr double keV
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:63
static constexpr double GeV
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
HepLorentzVector & boost(double, double, double)
void SetStatusChange(G4HadFinalStateStatus aS)
const G4ParticleDefinition * electron