| 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 // Created on 2016/04/08
27 //
28 // Authors: D. Sakata, S. Incerti
29 //
30 // This class perform transmission term of vol
31 // based on Quinn Model, see Phys. Rev. vol 12
32
33 #include "G4DNAQuinnPlasmonExcitationModel.hh"
34 #include "G4SystemOfUnits.hh"
35 #include "G4RandomDirection.hh"
36
37 //....oooOO0OOooo........oooOO0OOooo........oo
38
39 using namespace std;
40
41 //....oooOO0OOooo........oooOO0OOooo........oo
42
43 G4DNAQuinnPlasmonExcitationModel::G4DNAQuinnPl
44 (const G4ParticleDefin
45 const G4String& nam):
46 G4VEmModel(nam)
47 {
48 fpMaterialDensity = nullptr;
49 fLowEnergyLimit = 10 * eV;
50 fHighEnergyLimit = 1.0 * GeV;
51
52 for(int & i : nValenceElectron) i=0;
53
54 verboseLevel = 0;
55
56 if (verboseLevel > 0)
57 {
58 G4cout << "Quinn plasmon excitation model
59 }
60 fParticleChangeForGamma = nullptr;
61 statCode = false;
62 }
63
64 //....oooOO0OOooo........oooOO0OOooo........oo
65
66 G4DNAQuinnPlasmonExcitationModel::~G4DNAQuinnP
67 = default;
68
69 //....oooOO0OOooo........oooOO0OOooo........oo
70
71 void G4DNAQuinnPlasmonExcitationModel::Initial
72 (const G4ParticleDefin
73 const G4DataVector& /
74 {
75 for(int & i : nValenceElectron) i=0;
76
77 if (verboseLevel > 3)
78 {
79 G4cout <<
80 "Calling G4DNAQuinnPlasmonExcitatio
81 << G4endl;
82 }
83
84 if(particle == G4Electron::ElectronDefinitio
85 {
86 fLowEnergyLimit = 10 * eV;
87 fHighEnergyLimit = 1.0 * GeV;
88 }
89 else
90 {
91 G4Exception("G4DNAQuinnPlasmonExcitationMo
92 FatalException,"Not defined for other part
93 return;
94 }
95
96 // Get Number of valence electrons
97 G4ProductionCutsTable* theCoupleTable =
98 G4ProductionCutsTable::GetProductionCuts
99
100 auto numOfCouples = (G4int)theCoupleTable->
101
102 for(G4int i=0;i<numOfCouples;i++){
103
104 const G4MaterialCutsCouple* couple =
105 theCoupleTable->GetMaterialCutsCoupl
106
107 const G4Material* material = couple->GetMa
108
109 const G4ElementVector* theElementVector =m
110
111 std::size_t nelm = material->GetNumberOfEl
112 if (nelm==1) // Protection: only for singl
113 {
114 G4int z = G4lrint((*theElementVector)[0]
115 if(z<=100)
116 {
117 nValenceElectron[z] = GetNValenceElec
118 }
119 else
120 {
121 G4Exception("G4DNAQuinnPlasmonExcitati
122 FatalException,"The model is not app
123 }
124 }
125 //for(G4int j=0;j<nelm;j++){
126 // G4int z=G4lrint((*theElementVector)[
127 // if(z<=100){nValenceElectron[z] = Ge
128 //}
129 }
130
131 if( verboseLevel>0 )
132 {
133 G4cout << "Quinn plasmon excitation model
134 << "Energy range: "
135 << LowEnergyLimit() / eV << " eV - "
136 << HighEnergyLimit() / keV << " keV for "
137 << particle->GetParticleName()
138 << G4endl;
139 }
140
141 if (isInitialised){return;}
142 fParticleChangeForGamma = GetParticleChangeF
143 isInitialised = true;
144 }
145
146 //....oooOO0OOooo........oooOO0OOooo........oo
147
148 G4double G4DNAQuinnPlasmonExcitationModel::Cro
149 (const G4Material* ma
150 const G4ParticleDefi
151 G4double ekin,
152 G4double,
153 G4double)
154 {
155 if (verboseLevel > 3)
156 {
157 G4cout <<
158 "Calling CrossSectionPerVolume() of G
159 << G4endl;
160 }
161
162 // Protection: only for single element
163 if(material->GetNumberOfElements()>1) return
164 G4double z = material->GetZ();
165
166 // Protection: only for Gold
167 if (z!=79){return 0.;}
168
169
170 G4double sigma = 0;
171 G4double atomicNDensity = material->GetAtomi
172
173 if(atomicNDensity!= 0.0)
174 {
175 if (ekin >= fLowEnergyLimit && ekin < fHig
176 {
177 sigma = GetCrossSection(material,particl
178 }
179
180 if (verboseLevel > 2)
181 {
182 G4cout<<"_______________________________
183 G4cout<<"=== G4DNAQuinnPlasmonExcitation
184 G4cout<<"=== Kinetic energy (eV)=" << ek
185 <<particleDefinition->GetParticleN
186 G4cout<<"=== Cross section per atom for
187 <<sigma/cm/cm << G4endl;
188 G4cout<<"=== Cross section per atom for
189 <<sigma*atomicNDensity/(1./cm) <<
190 G4cout<<"=== G4DNAQuinnPlasmonExcitation
191 }
192 }
193
194 return sigma*atomicNDensity;
195 }
196
197 //....oooOO0OOooo........oooOO0OOooo........oo
198
199 void G4DNAQuinnPlasmonExcitationModel::SampleS
200 (std::vector<G4Dynami
201 const G4MaterialCuts
202 const G4DynamicParti
203 G4double,G4double)
204 {
205
206 if (verboseLevel > 3)
207 {
208 G4cout <<
209 "Calling SampleSecondaries() of G4D
210 << G4endl;
211 }
212
213 const G4Material *material = couple->GetMate
214
215 G4ParticleDefinition* particle = aDynamicPar
216
217 G4double k = aDynamicPar
218
219 if(particle == G4Electron::ElectronDefinitio
220 {
221 G4double e = 1.;
222 G4int z = material->GetZ();
223 G4int Nve = 0;
224
225 //TODO: have to be change to realistic!!
226 if(z<100) Nve = nValenceElectron[z];
227
228 G4double A = material->GetA()/g/mole
229 G4double Dens = material->GetDensity()/
230 G4double veDens = Dens*CLHEP::Avogadro*Nv
231
232 G4double omega_p = std::sqrt(veDens*std::p
233 (CLHEP::epsilon0/(1./cm)
234 /(CLHEP::c_squared/cm/cm
235
236 G4double excitationEnergy = CLHEP::hbar_
237 G4double newEnergy = k - excitati
238
239
240 if (newEnergy > 0)
241 {
242 fParticleChangeForGamma->
243 ProposeMomentumDirection(aDynamicP
244
245 fParticleChangeForGamma->ProposeLocalEne
246
247 if(!statCode)
248 {
249 fParticleChangeForGamma->SetPropose
250 }
251 else
252 {
253 fParticleChangeForGamma->SetPropose
254
255 }
256 }
257 }
258 }
259
260 //....oooOO0OOooo........oooOO0OOooo........oo
261
262 G4double G4DNAQuinnPlasmonExcitationModel::Get
263 (const G4Material* ma
264 const G4ParticleDefi
265 G4double kineticEner
266 {
267 G4double value=0;
268
269 if(particle == G4Electron::ElectronDefinitio
270 {
271 G4double e = 1.;
272 G4int z = material->GetZ();
273 G4int Nve = 0;
274 if(z<100) Nve = nValenceElectron[z];
275 G4double A = material->GetA()/g/mol
276 G4double Dens = material->GetDensity()
277 G4double veDens = Dens*CLHEP::Avogadro*N
278
279 G4double omega_p = std::sqrt(veDens*std::
280 /(CLHEP::epsilon0/(1./
281 (CLHEP::c_squared/cm/c
282
283 G4double fEnergy = std::pow(CLHEP::h_Plan
284 std::pow(3*veDens/CLHE
285 *(CLHEP::c_squared/cm/
286
287 G4double p0 = sqrt(2*CLHEP::electron
288 /(CLHEP::c_squared/cm/
289
290 G4double p = sqrt(2*CLHEP::electron
291 /(CLHEP::c_squared/cm/
292
293 G4double mfp = 2*CLHEP::Bohr_radius/c
294 /(CLHEP::hbar_Planck*o
295 (G4Log((std::pow(std::
296 +2*CLHEP::electron_mas
297 (CLHEP::c_squared/cm/c
298 *CLHEP::hbar_Planck,1.
299 /(p-std::pow(std::pow(
300 (CLHEP::c_squared/cm/c
301 *CLHEP::hbar_Planck,1.
302
303 G4double excitationEnergy = CLHEP::hbar
304
305 if((0<mfp)&&(0<veDens)&&(excitationEnergy
306 value = 1./(veDens*mfp);
307 }
308 }
309 return value*cm*cm;
310 }
311
312 //....oooOO0OOooo........oooOO0OOooo........oo
313
314 G4int G4DNAQuinnPlasmonExcitationModel::GetNVa
315 {
316
317 G4int Nve=0;
318
319 // Current limitation to gold
320 if (z!=79){return 0.;}
321
322 if (verboseLevel > 3)
323 {
324 G4cout <<
325 "Calling GetNValenceElectron() of G
326 << G4endl;
327 }
328
329 const char *datadir=nullptr;
330
331 if(datadir == nullptr)
332 {
333 datadir = G4FindDataDir("G4LEDATA");
334 if(datadir == nullptr)
335 {
336 G4Exception("G4DNAQuinnPlasmonExcitatio
337 ,"em0002",FatalException,
338 "Enviroment variable G4LEDA
339 return 0;
340 }
341 }
342
343 std::ostringstream targetfile;
344 targetfile.str("");
345 targetfile.clear(stringstream::goodbit);
346 targetfile << datadir <<"/dna/atomicstate_Z"
347 std::ifstream fin(targetfile.str().c_str());
348
349 if(!fin)
350 {
351 G4cout<< " Error : "<< targetfile.str() <<
352 G4Exception("G4DNAQuinnPlasmonExcitationMo
353 ,"em0003",FatalException,
354 "There is no target file");
355 return 0;
356 }
357
358 string buff0,buff1,buff2,buff3,buff4,buff5,b
359 fin >> buff0 >>buff1>>buff2>>buff3>>buff4>>b
360
361 while(true){
362 fin >> buff0 >>buff1>>buff2>>buff3>>buff4>
363 if(!fin.eof())
364 {
365 Nve = stoi(buff3);
366 }
367 else
368 {
369 break;
370 }
371 }
372 return Nve;
373 }
374
375