| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // Based on the work described in
27 // Rad Res 163, 98-111 (2005)
28 // D. Emfietzoglou, H. Nikjoo
29 //
30 // Authors of the class (2014):
31 // I. Kyriakou (kyriak@cc.uoi.gr)
32 // D. Emfietzoglou (demfietz@cc.uoi.gr)
33 // S. Incerti (incerti@cenbg.in2p3.fr)
34 //
35
36 #include "G4DNAEmfietzoglouExcitationModel.hh"
37 #include "G4SystemOfUnits.hh"
38 #include "G4DNAChemistryManager.hh"
39 #include "G4DNAMolecularMaterial.hh"
40
41 //....oooOO0OOooo........oooOO0OOooo........oo
42
43 using namespace std;
44
45 //....oooOO0OOooo........oooOO0OOooo........oo
46
47 G4DNAEmfietzoglouExcitationModel::G4DNAEmfietz
48
49 :G4VEmModel(nam)
50 {
51 fpMolWaterDensity = nullptr;
52
53 verboseLevel= 0;
54 // Verbosity scale:
55 // 0 = nothing
56 // 1 = warning for energy non-conservation
57 // 2 = details of energy budget
58 // 3 = calculation of cross sections, file
59 // 4 = entering in methods
60
61 if( verboseLevel>0 )
62 {
63 G4cout << "Emfietzoglou excitation model
64 }
65 fParticleChangeForGamma = nullptr;
66
67 SetLowEnergyLimit(8.*eV);
68 SetHighEnergyLimit(10.*keV);
69
70 // Selection of stationary mode
71 statCode = false;
72 }
73
74 //....oooOO0OOooo........oooOO0OOooo........oo
75
76 G4DNAEmfietzoglouExcitationModel::~G4DNAEmfiet
77 {
78 // Cross section
79
80 std::map< G4String,G4DNACrossSectionDataSe
81 for (pos = tableData.begin(); pos != table
82 {
83 G4DNACrossSectionDataSet* table = pos-
84 delete table;
85 }
86
87 }
88
89 //....oooOO0OOooo........oooOO0OOooo........oo
90
91 void G4DNAEmfietzoglouExcitationModel::Initial
92 cons
93 {
94
95 if (verboseLevel > 3)
96 G4cout << "Calling G4DNAEmfietzoglouEx
97
98 G4String fileElectron("dna/sigma_excitatio
99
100 G4ParticleDefinition* electronDef = G4Elec
101
102 G4String electron;
103
104 G4double scaleFactor = (1.e-22 / 3.343) *
105
106 // *** ELECTRON
107
108 electron = electronDef->GetParticleName();
109
110 tableFile[electron] = fileElectron;
111
112 // Cross section
113
114 auto tableE = new G4DNACrossSectionDataSe
115 tableE->LoadData(fileElectron);
116
117 tableData[electron] = tableE;
118
119 //
120
121 if( verboseLevel>0 )
122 {
123 G4cout << "Emfietzoglou excitation model
124 << "Energy range: "
125 << LowEnergyLimit() / eV << " eV
126 << HighEnergyLimit() / keV << " k
127 << particle->GetParticleName()
128 << G4endl;
129 }
130
131 // Initialize water density pointer
132 fpMolWaterDensity = G4DNAMolecularMaterial
133
134 if (isInitialised) return;
135 fParticleChangeForGamma = GetParticleChang
136 isInitialised = true;
137 }
138
139 //....oooOO0OOooo........oooOO0OOooo........oo
140
141 G4double G4DNAEmfietzoglouExcitationModel::Cro
142
143
144
145
146 {
147 if (verboseLevel > 3)
148 G4cout << "Calling CrossSectionPerVolu
149
150 if (particleDefinition != G4Electron::Elec
151
152 // Calculate total cross section for model
153
154 G4double sigma=0;
155
156 G4double waterDensity = (*fpMolWaterDensit
157
158 const G4String& particleName = particleDef
159
160 if (ekin >= LowEnergyLimit() && ekin <= Hi
161 {
162 std::map< G4String,G4DNACrossSectionData
163 pos = tableData.find(particleName);
164
165 if (pos != tableData.end())
166 {
167 G4DNACrossSectionDataSet* table = pos-
168 if (table != nullptr) sigma = table->F
169 }
170 else
171 {
172 G4Exception("G4DNAEmfietzoglouExcitati
173 FatalException,"Mo
174 }
175 }
176
177 if (verboseLevel > 2)
178 {
179 G4cout << "_____________________________
180 G4cout << "G4DNAEmfietzoglouExcitationMo
181 G4cout << "Kinetic energy(eV)=" << ekin/
182 G4cout << "Cross section per water molec
183 G4cout << "Cross section per water molec
184 //G4cout << " Cross section per water
185 ///sigma*material->GetAtomicNumDensityVe
186 G4cout << "G4DNAEmfietzoglouExcitationMo
187 }
188
189 return sigma*waterDensity;
190 }
191
192 //....oooOO0OOooo........oooOO0OOooo........oo
193
194 void G4DNAEmfietzoglouExcitationModel::SampleS
195
196
197
198
199 {
200
201 if (verboseLevel > 3)
202 G4cout << "Calling SampleSecondaries()
203
204 G4double k = aDynamicParticle->GetKineticE
205
206 const G4String& particleName = aDynamicPar
207
208 G4int level = RandomSelect(k,particleName)
209 G4double excitationEnergy = waterStructure
210 G4double newEnergy = k - excitationEnergy;
211
212 if (newEnergy > 0)
213 {
214 fParticleChangeForGamma->ProposeMoment
215
216 if (!statCode) fParticleChangeForGamma
217 else fParticleChangeForGamma->SetPropo
218
219 fParticleChangeForGamma->ProposeLocalE
220 }
221
222 const G4Track * theIncomingTrack = fPartic
223 G4DNAChemistryManager::Instance()->CreateW
224
225
226 }
227
228 //....oooOO0OOooo........oooOO0OOooo........oo
229
230 G4int G4DNAEmfietzoglouExcitationModel::Random
231 {
232
233 G4int level = 0;
234
235 std::map< G4String,G4DNACrossSectionDataSe
236 pos = tableData.find(particle);
237
238 if (pos != tableData.end())
239 {
240 G4DNACrossSectionDataSet* table = pos-
241
242 if (table != nullptr)
243 {
244 auto valuesBuffer = new G4double[
245 const auto n = (G4int)table->Numb
246 G4int i(n);
247 G4double value = 0.;
248
249 //Check reading of initial xs file
250 //G4cout << table->GetComponent(0)->
251 //G4cout << table->GetComponent(1)
252 //G4cout << table->GetComponent(2)
253 //G4cout << table->GetComponent(3)
254 //G4cout << table->GetComponent(4)
255 //G4cout << table->GetComponent(5)
256 //G4cout << table->GetComponent(6)
257 //abort();
258
259 while (i>0)
260 {
261 i--;
262 valuesBuffer[i] = table->GetCo
263 value += valuesBuffer[i];
264 }
265
266 value *= G4UniformRand();
267
268 i = n;
269
270 while (i > 0)
271 {
272 i--;
273
274 if (valuesBuffer[i] > value)
275 {
276 delete[] valuesBuffer;
277 return i;
278 }
279 value -= valuesBuffer[i];
280 }
281
282 delete[] valuesBuffer;
283
284 }
285 }
286 else
287 {
288 G4Exception("G4DNAEmfietzoglouExcitati
289 FatalException,"Model not
290 }
291 return level;
292 }
293