source: trunk/source/processes/electromagnetic/utils/src/G4EmSaturation.cc @ 1337

Last change on this file since 1337 was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales
  • Property svn:executable set to *
File size: 8.7 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// $Id: G4EmSaturation.cc,v 1.10 2009/09/25 09:16:40 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-04-beta-01 $
28//
29// -------------------------------------------------------------------
30//
31// GEANT4 Class file
32//
33//
34// File name:     G4EmSaturation
35//
36// Author:        Vladimir Ivanchenko
37//
38// Creation date: 18.02.2008
39//
40// Modifications:
41//
42// -------------------------------------------------------------
43
44//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
45//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
46
47#include "G4EmSaturation.hh"
48#include "G4LossTableManager.hh"
49#include "G4NistManager.hh"
50#include "G4Material.hh"
51#include "G4MaterialCutsCouple.hh"
52#include "G4Electron.hh"
53#include "G4Proton.hh"
54
55//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
56
57G4EmSaturation::G4EmSaturation()
58{
59  verbose = 1;
60  manager = 0;
61  curMaterial = 0;
62  curBirks = 0.0;
63  curRatio = 1.0;
64  curChargeSq = 1.0;
65  nMaterials = 0;
66  electron = 0;
67  Initialise(); 
68}
69
70//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
71
72G4EmSaturation::~G4EmSaturation()
73{}
74
75//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
76
77G4double G4EmSaturation::VisibleEnergyDeposition(
78                                      const G4ParticleDefinition* p, 
79                                      const G4MaterialCutsCouple* couple, 
80                                      G4double length,
81                                      G4double edep,
82                                      G4double niel)
83{
84  if(edep <= 0.0) return 0.0;
85
86  G4double evis = edep;
87  G4double bfactor = FindBirksCoefficient(couple->GetMaterial());
88
89  if(bfactor > 0.0) { 
90
91    G4int pdgCode = p->GetPDGEncoding();
92    // atomic relaxations for gamma incident
93    if(22 == pdgCode) {
94      evis /= (1.0 + bfactor*edep/manager->GetRange(electron,edep,couple));
95
96      // energy loss
97    } else {
98
99      // protections
100      G4double nloss = niel;
101      if(nloss < 0.0) nloss = 0.0;
102      G4double eloss = edep - nloss;
103
104      // neutrons
105      if(2112 == pdgCode || eloss < 0.0 || length <= 0.0) {
106        nloss = edep;
107        eloss = 0.0;
108      }
109
110      // continues energy loss
111      if(eloss > 0.0) eloss /= (1.0 + bfactor*eloss/length);
112 
113      // non-ionizing energy loss
114      if(nloss > 0.0) {
115        if(!proton) {proton = G4Proton::Proton();}
116        G4double escaled = nloss*curRatio;
117        G4double s = manager->GetRange(proton,escaled,couple)/curChargeSq; 
118        nloss /= (1.0 + bfactor*nloss/s);
119      }
120
121      evis = eloss + nloss;
122    }
123  }
124 
125  return evis;
126}
127
128//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
129
130G4double G4EmSaturation::FindG4BirksCoefficient(const G4Material* mat)
131{
132  G4String name = mat->GetName();
133  // is this material in the vector?
134 
135  for(G4int j=0; j<nG4Birks; j++) {
136    if(name == g4MatNames[j]) {
137      if(verbose > 0) 
138        G4cout << "### G4EmSaturation::FindG4BirksCoefficient for "
139               << name << " is " << g4MatData[j]*MeV/mm << " mm/MeV "
140               << G4endl;
141      return g4MatData[j];
142    }
143  }
144  return FindBirksCoefficient(mat);
145}
146
147//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
148
149G4double G4EmSaturation::FindBirksCoefficient(const G4Material* mat)
150{
151  // electron should exist in any case
152  if(!manager) {
153    manager = G4LossTableManager::Instance();
154    nist    = G4NistManager::Instance();
155    electron= G4Electron::Electron();
156    proton  = 0;
157  }
158
159  if(mat == curMaterial) return curBirks;
160
161  curMaterial = mat;
162  curBirks = 0.0;
163  curRatio = 1.0;
164  curChargeSq = 1.0;
165
166  // seach in the run-time list
167  for(G4int i=0; i<nMaterials; i++) {
168    if(mat == matPointers[i]) {
169      curBirks = mat->GetIonisation()->GetBirksConstant();
170      curRatio = massFactors[i];
171      curChargeSq = effCharges[i];
172      return curBirks;
173    }
174  }
175
176  G4String name = mat->GetName();
177  curBirks = mat->GetIonisation()->GetBirksConstant();
178
179  // material has no Birks coeffitient defined
180  // seach in the Geant4 list
181  if(curBirks == 0.0) {
182    for(G4int j=0; j<nG4Birks; j++) {
183      if(name == g4MatNames[j]) {
184        mat->GetIonisation()->SetBirksConstant(g4MatData[j]);
185        curBirks = g4MatData[j];
186        break;
187      }
188    }
189  }
190
191  if(curBirks == 0.0 && verbose > 0) {
192      G4cout << "### G4EmSaturation::FindBirksCoefficient fails "
193        " for material " << name << G4endl;
194  }
195
196  // compute mean mass ratio
197  curRatio = 0.0;
198  curChargeSq = 0.0;
199  G4double norm = 0.0;
200  const G4ElementVector* theElementVector = mat->GetElementVector();
201  const G4double* theAtomNumDensityVector = mat->GetVecNbOfAtomsPerVolume();
202  size_t nelm = mat->GetNumberOfElements();
203  for (size_t i=0; i<nelm; i++) {
204    const G4Element* elm = (*theElementVector)[i];
205    G4double Z = elm->GetZ();
206    G4double w = Z*Z*theAtomNumDensityVector[i];
207    curRatio += w/nist->GetAtomicMassAmu(G4int(Z));
208    curChargeSq = Z*Z*w;
209    norm += w;
210  }
211  curRatio *= proton_mass_c2/norm;
212  curChargeSq /= norm;
213
214  // store results
215  matPointers.push_back(mat);
216  matNames.push_back(name);
217  massFactors.push_back(curRatio);
218  effCharges.push_back(curChargeSq);
219  nMaterials++;
220  if(curBirks > 0.0 && verbose > 0) {
221    G4cout << "### G4EmSaturation::FindBirksCoefficient Birks coefficient for "
222           << name << "  " << curBirks*MeV/mm << " mm/MeV" << G4endl;
223  }
224  return curBirks;
225}
226
227//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
228
229void G4EmSaturation::DumpBirksCoefficients()
230{
231  if(nMaterials > 0) {
232    G4cout << "### Birks coeffitients used in run time" << G4endl;
233    for(G4int i=0; i<nMaterials; i++) {
234      G4double br = matPointers[i]->GetIonisation()->GetBirksConstant();
235      G4cout << "   " << matNames[i] << "     " 
236             << br*MeV/mm << " mm/MeV" << "     "
237             << br*matPointers[i]->GetDensity()*MeV*cm2/g
238             << " g/cm^2/MeV" 
239             << G4endl;
240    }
241  }
242}
243
244//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
245
246void G4EmSaturation::DumpG4BirksCoefficients()
247{
248  if(nG4Birks > 0) {
249    G4cout << "### Birks coeffitients for Geant4 materials" << G4endl;
250    for(G4int i=0; i<nG4Birks; i++) {
251      G4cout << "   " << g4MatNames[i] << "   " 
252             << g4MatData[i]*MeV/mm << " mm/MeV" << G4endl;
253    }
254  }
255}
256
257//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
258
259void G4EmSaturation::Initialise()
260{
261  // M.Hirschberg et al., IEEE Trans. Nuc. Sci. 39 (1992) 511
262  // SCSN-38 kB = 0.00842 g/cm^2/MeV; rho = 1.06 g/cm^3
263  g4MatNames.push_back("G4_POLYSTYRENE");
264  g4MatData.push_back(0.07943*mm/MeV);
265
266  // C.Fabjan (private communication)
267  // kB = 0.006 g/cm^2/MeV; rho = 7.13 g/cm^3
268  g4MatNames.push_back("G4_BGO");
269  g4MatData.push_back(0.008415*mm/MeV);
270
271  // A.Ribon analysis of publications
272  // Scallettar et al., Phys. Rev. A25 (1982) 2419.
273  // NIM A 523 (2004) 275.
274  // kB = 0.022 g/cm^2/MeV; rho = 1.396 g/cm^3;
275  // ATLAS Efield = 10 kV/cm provide the strongest effect
276  g4MatNames.push_back("G4_lAr");
277  g4MatData.push_back(0.1576*mm/MeV);
278
279  //G4_BARIUM_FLUORIDE
280  //G4_CESIUM_IODIDE
281  //G4_GEL_PHOTO_EMULSION
282  //G4_PHOTO_EMULSION
283  //G4_PLASTIC_SC_VINYLTOLUENE
284  //G4_SODIUM_IODIDE
285  //G4_STILBENE
286  //G4_lAr
287  //G4_PbWO4
288  //G4_Lucite
289
290  nG4Birks = g4MatData.size();
291}
292
293//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
Note: See TracBrowser for help on using the repository browser.