source: trunk/source/particles/management/src/G4NucleiProperties.cc@ 1192

Last change on this file since 1192 was 992, checked in by garnier, 17 years ago

fichiers oublies
File size: 8.5 KB
Line 
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//
27// $Id: G4NucleiProperties.cc,v 1.18 2008/11/06 13:17:36 kurasige Exp $
28// GEANT4 tag $Name: geant4-09-02-ref-02 $
29//
30//
31// ------------------------------------------------------------
32// GEANT 4 class header file
33//
34// ------------------------------------------------------------
35//
36// Hadronic Process: Nuclear De-excitations
37// by V. Lara (Oct 1998)
38// Migrate into particles category by H.Kurashige (17 Nov. 98)
39// Added Shell-Pairing corrections to the Cameron mass
40// excess formula by V.Lara (9 May 99)
41
42
43#include "G4NucleiProperties.hh"
44
45G4bool G4NucleiProperties::isIntialized = false;
46
47G4double G4NucleiProperties::mass_proton = -1.;
48G4double G4NucleiProperties::mass_neutron = -1.;
49G4double G4NucleiProperties::mass_deuteron = -1.;
50G4double G4NucleiProperties::mass_triton = -1.;
51G4double G4NucleiProperties::mass_alpha = -1.;
52G4double G4NucleiProperties::mass_He3 = -1.;
53G4double G4NucleiProperties::electronMass[MaxZ];
54
55G4double G4NucleiProperties::AtomicMass(G4double A, G4double Z)
56{
57 const G4double hydrogen_mass_excess = G4NucleiPropertiesTable::GetMassExcess(1,1);
58 const G4double neutron_mass_excess = G4NucleiPropertiesTable::GetMassExcess(0,1);
59
60 G4double mass =
61 (A-Z)*neutron_mass_excess + Z*hydrogen_mass_excess - BindingEnergy(A,Z) + A*amu_c2;
62
63 return mass;
64}
65
66G4double G4NucleiProperties::BindingEnergy(G4double A, G4double Z)
67{
68 //
69 // Weitzsaecker's Mass formula
70 //
71 G4int Npairing = G4int(A-Z)%2; // pairing
72 G4int Zpairing = G4int(Z)%2;
73 G4double binding =
74 - 15.67*A // nuclear volume
75 + 17.23*std::pow(A,2./3.) // surface energy
76 + 93.15*((A/2.-Z)*(A/2.-Z))/A // asymmetry
77 + 0.6984523*Z*Z*std::pow(A,-1./3.); // coulomb
78 if( Npairing == Zpairing ) binding += (Npairing+Zpairing-1) * 12.0 / std::sqrt(A); // pairing
79
80 return -binding*MeV;
81}
82
83
84G4double G4NucleiProperties::GetNuclearMass(const G4double A, const G4double Z)
85{
86 if (!isIntialized) {
87 isIntialized = true;
88 G4ParticleDefinition * nucleus = 0;
89 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("proton"); // proton
90 if (nucleus!=0) mass_proton = nucleus->GetPDGMass();
91 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("neutron"); // neutron
92 if (nucleus!=0) mass_neutron = nucleus->GetPDGMass();
93 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("deuteron"); // deuteron
94 if (nucleus!=0) mass_deuteron = nucleus->GetPDGMass();
95 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("triton"); // triton
96 if (nucleus!=0) mass_triton = nucleus->GetPDGMass();
97 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("alpha"); // alpha
98 if (nucleus!=0) mass_alpha = nucleus->GetPDGMass();
99 nucleus = G4ParticleTable::GetParticleTable()->FindParticle("He3"); // He3
100 if (nucleus!=0) mass_He3 = nucleus->GetPDGMass();
101
102 for (int iz=1; iz<MaxZ; iz+=1){
103 electronMass[iz] = iz*electron_mass_c2
104 - 1.433e-5*MeV*std::pow(G4double(iz),2.39); ;
105 }
106 }
107
108 if (A < 1 || Z < 0 || Z > A) {
109#ifdef G4VERBOSE
110 if (G4ParticleTable::GetParticleTable()->GetVerboseLevel()>0) {
111 G4cout << "G4NucleiProperties::GetNuclearMass: Wrong values for A = " << A
112 << " and Z = " << Z << G4endl;
113 }
114#endif
115 return 0.0;
116 } else {
117 G4double mass= -1.;
118 if ( (Z<=2) ) {
119 if ( (Z==1)&&(A==1) ) {
120 mass = mass_proton;
121 } else if ( (Z==0)&&(A==1) ) {
122 mass = mass_neutron;
123 } else if ( (Z==1)&&(A==2) ) {
124 mass = mass_deuteron;
125 } else if ( (Z==1)&&(A==3) ) {
126 mass = mass_triton;
127 } else if ( (Z==2)&&(A==4) ) {
128 mass = mass_alpha;
129 } else if ( (Z==2)&&(A==3) ) {
130 mass = mass_He3;
131 }
132 }
133 if (mass < 0.) {
134 if (Z >= MaxZ) {
135 mass = GetAtomicMass(A,Z) - Z*electron_mass_c2 + 1.433e-5*MeV*std::pow(Z,2.39);
136 } else {
137 mass = GetAtomicMass(A,Z) - electronMass[G4int(Z)];
138 }
139 }
140 if (mass < 0.) mass = 0.0;
141 return mass;
142 }
143}
144
145G4bool G4NucleiProperties::IsInStableTable(const G4double A, const G4double Z)
146{
147 if (Z < 0 || Z > A) {
148#ifdef G4VERBOSE
149 if (G4ParticleTable::GetParticleTable()->GetVerboseLevel()>0) {
150 G4cout << "G4NucleiProperties::IsInStableTable: Wrong values for A = "
151 << A << " and Z = " << Z << G4endl;
152 }
153#endif
154 return false;
155
156 } else {
157 G4int iA = G4int(A);
158 G4int iZ = G4int(Z);
159 return G4NucleiPropertiesTable::IsInTable(iZ,iA);
160 }
161}
162
163G4double G4NucleiProperties::GetMassExcess(const G4double A, const G4double Z)
164{
165 G4int iA = G4int(A);
166 G4int iZ = G4int(Z);
167 return GetMassExcess(iA,iZ);
168}
169
170G4double G4NucleiProperties::GetMassExcess(const G4int A, const G4int Z)
171{
172 if (A < 1 || Z < 0 || Z > A) {
173#ifdef G4VERBOSE
174 if (G4ParticleTable::GetParticleTable()->GetVerboseLevel()>0) {
175 G4cout << "G4NucleiProperties::GetMassExccess: Wrong values for A = "
176 << A << " and Z = " << Z << G4endl;
177 }
178#endif
179 return 0.0;
180
181 } else {
182
183 if (G4NucleiPropertiesTable::IsInTable(Z,A)){
184 return G4NucleiPropertiesTable::GetMassExcess(Z,A);
185 } else if (G4NucleiPropertiesTheoreticalTable::IsInTable(Z,A)){
186 return G4NucleiPropertiesTheoreticalTable::GetMassExcess(Z,A);
187 } else {
188 return MassExcess(A,Z);
189 }
190 }
191
192}
193
194
195G4double G4NucleiProperties::GetAtomicMass(const G4double A, const G4double Z)
196{
197 if (Z < 0 || Z > A) {
198#ifdef G4VERBOSE
199 if (G4ParticleTable::GetParticleTable()->GetVerboseLevel()>0) {
200 G4cout << "G4NucleiProperties::GetAtomicMass: Wrong values for A = "
201 << A << " and Z = " << Z << G4endl;
202 }
203#endif
204 return 0.0;
205
206 } else if (std::abs(A - G4int(A)) > 1.e-10) {
207 return AtomicMass(A,Z);
208
209 } else {
210 G4int iA = G4int(A);
211 G4int iZ = G4int(Z);
212 if (G4NucleiPropertiesTable::IsInTable(iZ,iA)) {
213 return G4NucleiPropertiesTable::GetAtomicMass(iZ,iA);
214 } else if (G4NucleiPropertiesTheoreticalTable::IsInTable(iZ,iA)){
215 return G4NucleiPropertiesTheoreticalTable::GetAtomicMass(iZ,iA);
216 } else {
217 return AtomicMass(A,Z);
218 }
219 }
220}
221
222G4double G4NucleiProperties::GetBindingEnergy(const G4double A, const G4double Z)
223{
224 G4int iA = G4int(A);
225 G4int iZ = G4int(Z);
226 return GetBindingEnergy(iA,iZ);
227}
228
229G4double G4NucleiProperties::GetBindingEnergy(const G4int A, const G4int Z)
230{
231 if (A < 1 || Z < 0 || Z > A) {
232#ifdef G4VERBOSE
233 if (G4ParticleTable::GetParticleTable()->GetVerboseLevel()>0) {
234 G4cout << "G4NucleiProperties::GetMassExccess: Wrong values for A = "
235 << A << " and Z = " << Z << G4endl;
236 }
237#endif
238 return 0.0;
239
240 } else {
241 if (G4NucleiPropertiesTable::IsInTable(Z,A)) {
242 return G4NucleiPropertiesTable::GetBindingEnergy(Z,A);
243 } else if (G4NucleiPropertiesTheoreticalTable::IsInTable(Z,A)) {
244 return G4NucleiPropertiesTheoreticalTable::GetBindingEnergy(Z,A);
245 }else {
246 return BindingEnergy(A,Z);
247 }
248
249 }
250}
251
252
253G4double G4NucleiProperties::MassExcess(G4double A, G4double Z)
254{
255 return GetAtomicMass(A,Z) - A*amu_c2;
256}
257
Note: See TracBrowser for help on using the repository browser.