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

Last change on this file since 1117 was 992, checked in by garnier, 15 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.