source: trunk/source/processes/electromagnetic/utils/src/G4VEmModel.cc@ 1330

Last change on this file since 1330 was 1315, checked in by garnier, 15 years ago

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08

File size: 9.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// $Id: G4VEmModel.cc,v 1.35 2010/05/26 18:06:25 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
28//
29// -------------------------------------------------------------------
30//
31// GEANT4 Class file
32//
33//
34// File name: G4VEmModel
35//
36// Author: Vladimir Ivanchenko
37//
38// Creation date: 25.07.2005
39//
40// Modifications:
41// 25.10.2005 Set default highLimit=100.TeV (V.Ivanchenko)
42// 06.02.2006 add method ComputeMeanFreePath() (mma)
43// 16.02.2009 Move implementations of virtual methods to source (VI)
44//
45//
46// Class Description:
47//
48// Abstract interface to energy loss models
49
50// -------------------------------------------------------------------
51//
52
53#include "G4VEmModel.hh"
54#include "G4LossTableManager.hh"
55#include "G4ProductionCutsTable.hh"
56#include "G4ParticleChangeForLoss.hh"
57#include "G4ParticleChangeForGamma.hh"
58
59//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
60//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
61
62G4VEmModel::G4VEmModel(const G4String& nam):
63 fluc(0), name(nam), lowLimit(0.1*keV), highLimit(100.0*TeV),
64 eMinActive(0.0),eMaxActive(DBL_MAX),
65 polarAngleLimit(0.0),secondaryThreshold(DBL_MAX),theLPMflag(false),
66 pParticleChange(0),nuclearStopping(false),
67 currentCouple(0),currentElement(0),
68 nsec(5),flagDeexcitation(false)
69{
70 xsec.resize(nsec);
71 nSelectors = 0;
72 G4LossTableManager::Instance()->Register(this);
73}
74
75//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
76
77G4VEmModel::~G4VEmModel()
78{
79 G4LossTableManager::Instance()->DeRegister(this);
80 G4int n = elmSelectors.size();
81 if(n > 0) {
82 for(G4int i=0; i<n; ++i) {
83 delete elmSelectors[i];
84 }
85 }
86}
87
88//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
89
90G4ParticleChangeForLoss* G4VEmModel::GetParticleChangeForLoss()
91{
92 G4ParticleChangeForLoss* p = 0;
93 if (pParticleChange) {
94 p = static_cast<G4ParticleChangeForLoss*>(pParticleChange);
95 } else {
96 p = new G4ParticleChangeForLoss();
97 pParticleChange = p;
98 }
99 return p;
100}
101
102//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
103
104G4ParticleChangeForGamma* G4VEmModel::GetParticleChangeForGamma()
105{
106 G4ParticleChangeForGamma* p = 0;
107 if (pParticleChange) {
108 p = static_cast<G4ParticleChangeForGamma*>(pParticleChange);
109 } else {
110 p = new G4ParticleChangeForGamma();
111 pParticleChange = p;
112 }
113 return p;
114}
115
116//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
117
118void G4VEmModel::InitialiseElementSelectors(const G4ParticleDefinition* p,
119 const G4DataVector& cuts)
120{
121 // initialise before run
122 flagDeexcitation = false;
123 G4LossTableManager* man = G4LossTableManager::Instance();
124 G4bool spline = man->SplineFlag();
125
126 // two times less bins because probability functon is normalized
127 // so correspondingly is more smooth
128 G4int nbins = (man->GetNumberOfBinsPerDecade()/3)*
129 G4int(std::log10(highLimit/lowLimit) + 0.5);
130 if(nbins < 5) { nbins = 5; }
131
132 G4ProductionCutsTable* theCoupleTable=
133 G4ProductionCutsTable::GetProductionCutsTable();
134 G4int numOfCouples = theCoupleTable->GetTableSize();
135
136 // prepare vector
137 if(numOfCouples > nSelectors) {
138 elmSelectors.reserve(numOfCouples);
139 for(G4int i=nSelectors; i<numOfCouples; ++i) { elmSelectors.push_back(0); }
140 nSelectors = numOfCouples;
141 }
142
143 // initialise vector
144 for(G4int i=0; i<numOfCouples; ++i) {
145 currentCouple = theCoupleTable->GetMaterialCutsCouple(i);
146 const G4Material* material = currentCouple->GetMaterial();
147 G4int idx = currentCouple->GetIndex();
148
149 // selector already exist check if should be deleted
150 G4bool create = true;
151 if(elmSelectors[i]) {
152 if(material == elmSelectors[i]->GetMaterial()) { create = false; }
153 else { delete elmSelectors[i]; }
154 }
155 if(create) {
156 elmSelectors[i] = new G4EmElementSelector(this,material,nbins,
157 lowLimit,highLimit,spline);
158 }
159 elmSelectors[i]->Initialise(p, cuts[idx]);
160 //elmSelectors[i]->Dump(p);
161 }
162}
163
164//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
165
166G4double G4VEmModel::ComputeDEDXPerVolume(const G4Material*,
167 const G4ParticleDefinition*,
168 G4double,G4double)
169{
170 return 0.0;
171}
172
173//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
174
175G4double G4VEmModel::CrossSectionPerVolume(const G4Material* material,
176 const G4ParticleDefinition* p,
177 G4double ekin,
178 G4double emin,
179 G4double emax)
180{
181 SetupForMaterial(p, material, ekin);
182 G4double cross = 0.0;
183 const G4ElementVector* theElementVector = material->GetElementVector();
184 const G4double* theAtomNumDensityVector = material->GetVecNbOfAtomsPerVolume();
185 G4int nelm = material->GetNumberOfElements();
186 if(nelm > nsec) {
187 xsec.resize(nelm);
188 nsec = nelm;
189 }
190 for (G4int i=0; i<nelm; i++) {
191 cross += theAtomNumDensityVector[i]*
192 ComputeCrossSectionPerAtom(p,(*theElementVector)[i],ekin,emin,emax);
193 xsec[i] = cross;
194 }
195 return cross;
196}
197
198//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
199
200const G4Element* G4VEmModel::SelectRandomAtom(const G4Material* material,
201 const G4ParticleDefinition* pd,
202 G4double kinEnergy,
203 G4double tcut,
204 G4double tmax)
205{
206 const G4ElementVector* theElementVector = material->GetElementVector();
207 G4int n = material->GetNumberOfElements() - 1;
208 currentElement = (*theElementVector)[n];
209 if (n > 0) {
210 G4double x = G4UniformRand()*
211 G4VEmModel::CrossSectionPerVolume(material,pd,kinEnergy,tcut,tmax);
212 for(G4int i=0; i<n; ++i) {
213 if (x <= xsec[i]) {
214 currentElement = (*theElementVector)[i];
215 break;
216 }
217 }
218 }
219 return currentElement;
220}
221
222//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
223
224G4double G4VEmModel::ComputeCrossSectionPerAtom(const G4ParticleDefinition*,
225 G4double, G4double, G4double,
226 G4double, G4double)
227{
228 return 0.0;
229}
230
231//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
232
233void G4VEmModel::DefineForRegion(const G4Region*)
234{}
235
236//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
237
238G4double G4VEmModel::MinEnergyCut(const G4ParticleDefinition*,
239 const G4MaterialCutsCouple*)
240{
241 return 0.0;
242}
243
244//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
245
246G4double G4VEmModel::ChargeSquareRatio(const G4Track& track)
247{
248 return GetChargeSquareRatio(track.GetDefinition(),
249 track.GetMaterial(), track.GetKineticEnergy());
250}
251
252//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
253
254G4double G4VEmModel::GetChargeSquareRatio(const G4ParticleDefinition* p,
255 const G4Material*, G4double)
256{
257 G4double q = p->GetPDGCharge()/CLHEP::eplus;
258 return q*q;
259}
260
261//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
262
263G4double G4VEmModel::GetParticleCharge(const G4ParticleDefinition* p,
264 const G4Material*, G4double)
265{
266 return p->GetPDGCharge();
267}
268
269//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
270
271void G4VEmModel::CorrectionsAlongStep(const G4MaterialCutsCouple*,
272 const G4DynamicParticle*,
273 G4double&,G4double&,G4double)
274{}
275
276//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
277
278void G4VEmModel::SampleDeexcitationAlongStep(const G4Material*,
279 const G4Track&,
280 G4double& )
281{}
282
283//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
284
285G4double G4VEmModel::MaxSecondaryEnergy(const G4ParticleDefinition*,
286 G4double kineticEnergy)
287{
288 return kineticEnergy;
289}
290
291//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
292
293void G4VEmModel::SetupForMaterial(const G4ParticleDefinition*,
294 const G4Material*, G4double)
295{}
296
297//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
Note: See TracBrowser for help on using the repository browser.