Geant4  10.02
G4DNABornIonisationModel1.cc
Go to the documentation of this file.
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: G4DNABornIonisationModel1.cc 87631 2014-12-14 12:42:05Z matkara $
27 //
28 
29 #include "G4DNABornIonisationModel1.hh"
30 #include "G4PhysicalConstants.hh"
31 #include "G4SystemOfUnits.hh"
32 #include "G4UAtomicDeexcitation.hh"
33 #include "G4LossTableManager.hh"
34 #include "G4DNAChemistryManager.hh"
35 #include "G4DNAMolecularMaterial.hh"
36 #include "G4DNABornAngle.hh"
37 #include "G4DeltaAngle.hh"
38 
39 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
40 
41 using namespace std;
42 
43 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
44 
45 G4DNABornIonisationModel1::G4DNABornIonisationModel1(const G4ParticleDefinition*,
46  const G4String& nam) :
47  G4VEmModel(nam), isInitialised(false)
48 {
49  verboseLevel = 0;
50  // Verbosity scale:
51  // 0 = nothing
52  // 1 = warning for energy non-conservation
53  // 2 = details of energy budget
54  // 3 = calculation of cross sections, file openings, sampling of atoms
55  // 4 = entering in methods
56 
57  if (verboseLevel > 0)
58  {
59  G4cout << "Born ionisation model is constructed " << G4endl;
60  }
61 
62  //Mark this model as "applicable" for atomic deexcitation
63  SetDeexcitationFlag(true);
64  fAtomDeexcitation = 0;
65  fParticleChangeForGamma = 0;
66  fpMolWaterDensity = 0;
67 
68  // define default angular generator
69  SetAngularDistribution(new G4DNABornAngle());
70 
71  // Selection of computation method
72 
73  fasterCode = false;
74 
75  // Selection of stationary mode
76 
77  statCode = false;
78 
79  // Selection of SP scaling
80 
81  spScaling = true;
82 }
83 
84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
85 
86 G4DNABornIonisationModel1::~G4DNABornIonisationModel1()
87 {
88  // Cross section
89 
90  std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator pos;
91  for (pos = tableData.begin(); pos != tableData.end(); ++pos)
92  {
93  G4DNACrossSectionDataSet* table = pos->second;
94  delete table;
95  }
96 
97  // Final state
98 
99  eVecm.clear();
100  pVecm.clear();
101 }
102 
103 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
104 
105 void G4DNABornIonisationModel1::Initialise(const G4ParticleDefinition* particle,
106  const G4DataVector& /*cuts*/)
107 {
108 
109  if (verboseLevel > 3)
110  {
111  G4cout << "Calling G4DNABornIonisationModel1::Initialise()" << G4endl;
112  }
113 
114  // Energy limits
115 
116  G4String fileElectron("dna/sigma_ionisation_e_born");
117  G4String fileProton("dna/sigma_ionisation_p_born");
118 
119  G4ParticleDefinition* electronDef = G4Electron::ElectronDefinition();
120  G4ParticleDefinition* protonDef = G4Proton::ProtonDefinition();
121 
122  G4String electron;
123  G4String proton;
124 
125  G4double scaleFactor = (1.e-22 / 3.343) * m*m;
126 
127  char *path = getenv("G4LEDATA");
128 
129  // *** ELECTRON
130 
131  electron = electronDef->GetParticleName();
132 
133  tableFile[electron] = fileElectron;
134 
135  lowEnergyLimit[electron] = 11. * eV;
136  highEnergyLimit[electron] = 1. * MeV;
137 
138  // Cross section
139 
140  G4DNACrossSectionDataSet* tableE = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, eV,scaleFactor );
141  tableE->LoadData(fileElectron);
142 
143  tableData[electron] = tableE;
144 
145  // Final state
146 
147  std::ostringstream eFullFileName;
148 
149  if (fasterCode) eFullFileName << path << "/dna/sigmadiff_cumulated_ionisation_e_born_hp.dat";
150  if (!fasterCode) eFullFileName << path << "/dna/sigmadiff_ionisation_e_born.dat";
151 
152  std::ifstream eDiffCrossSection(eFullFileName.str().c_str());
153 
154  if (!eDiffCrossSection)
155  {
156  if (fasterCode) G4Exception("G4DNABornIonisationModel1::Initialise","em0003",
157  FatalException,"Missing data file:/dna/sigmadiff_cumulated_ionisation_e_born_hp.dat");
158 
159  if (!fasterCode) G4Exception("G4DNABornIonisationModel1::Initialise","em0003",
160  FatalException,"Missing data file:/dna/sigmadiff_ionisation_e_born.dat");
161  }
162 
163  //
164 
165  // Clear the arrays for re-initialization case (MT mode)
166  // March 25th, 2014 - Vaclav Stepan, Sebastien Incerti
167 
168  eTdummyVec.clear();
169  pTdummyVec.clear();
170 
171  eVecm.clear();
172  pVecm.clear();
173 
174  for (int j=0; j<5; j++)
175  {
176  eProbaShellMap[j].clear();
177  pProbaShellMap[j].clear();
178 
179  eDiffCrossSectionData[j].clear();
180  pDiffCrossSectionData[j].clear();
181 
182  eNrjTransfData[j].clear();
183  pNrjTransfData[j].clear();
184  }
185  //
186 
187  eTdummyVec.push_back(0.);
188  while(!eDiffCrossSection.eof())
189  {
190  double tDummy;
191  double eDummy;
192  eDiffCrossSection>>tDummy>>eDummy;
193  if (tDummy != eTdummyVec.back()) eTdummyVec.push_back(tDummy);
194 
195  double tmp;
196  for (int j=0; j<5; j++)
197  {
198  eDiffCrossSection>> tmp;
199 
200  eDiffCrossSectionData[j][tDummy][eDummy] = tmp;
201 
202  if (fasterCode)
203  {
204  eNrjTransfData[j][tDummy][eDiffCrossSectionData[j][tDummy][eDummy]]=eDummy;
205  eProbaShellMap[j][tDummy].push_back(eDiffCrossSectionData[j][tDummy][eDummy]);
206  }
207 
208  // SI - only if eof is not reached
209  if (!eDiffCrossSection.eof() && !fasterCode) eDiffCrossSectionData[j][tDummy][eDummy]*=scaleFactor;
210 
211  if (!fasterCode) eVecm[tDummy].push_back(eDummy);
212 
213  }
214  }
215 
216  // *** PROTON
217 
218  proton = protonDef->GetParticleName();
219 
220  tableFile[proton] = fileProton;
221 
222  lowEnergyLimit[proton] = 500. * keV;
223  highEnergyLimit[proton] = 100. * MeV;
224 
225  // Cross section
226 
227  G4DNACrossSectionDataSet* tableP = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, eV,scaleFactor );
228  tableP->LoadData(fileProton);
229 
230  tableData[proton] = tableP;
231 
232  // Final state
233 
234  std::ostringstream pFullFileName;
235 
236  if (fasterCode) pFullFileName << path << "/dna/sigmadiff_cumulated_ionisation_p_born_hp.dat";
237 
238  if (!fasterCode) pFullFileName << path << "/dna/sigmadiff_ionisation_p_born.dat";
239 
240  std::ifstream pDiffCrossSection(pFullFileName.str().c_str());
241 
242  if (!pDiffCrossSection)
243  {
244  if (fasterCode) G4Exception("G4DNABornIonisationModel1::Initialise","em0003",
245  FatalException,"Missing data file:/dna/sigmadiff_cumulated_ionisation_p_born_hp.dat");
246 
247  if (!fasterCode) G4Exception("G4DNABornIonisationModel1::Initialise","em0003",
248  FatalException,"Missing data file:/dna/sigmadiff_ionisation_p_born.dat");
249  }
250 
251  pTdummyVec.push_back(0.);
252  while(!pDiffCrossSection.eof())
253  {
254  double tDummy;
255  double eDummy;
256  pDiffCrossSection>>tDummy>>eDummy;
257  if (tDummy != pTdummyVec.back()) pTdummyVec.push_back(tDummy);
258  for (int j=0; j<5; j++)
259  {
260  pDiffCrossSection>>pDiffCrossSectionData[j][tDummy][eDummy];
261 
262  if (fasterCode)
263  {
264  pNrjTransfData[j][tDummy][pDiffCrossSectionData[j][tDummy][eDummy]]=eDummy;
265  pProbaShellMap[j][tDummy].push_back(pDiffCrossSectionData[j][tDummy][eDummy]);
266  }
267 
268  // SI - only if eof is not reached !
269  if (!pDiffCrossSection.eof() && !fasterCode) pDiffCrossSectionData[j][tDummy][eDummy]*=scaleFactor;
270 
271  if (!fasterCode) pVecm[tDummy].push_back(eDummy);
272  }
273  }
274 
275  //
276 
277  if (particle==electronDef)
278  {
279  SetLowEnergyLimit(lowEnergyLimit[electron]);
280  SetHighEnergyLimit(highEnergyLimit[electron]);
281  }
282 
283  if (particle==protonDef)
284  {
285  SetLowEnergyLimit(lowEnergyLimit[proton]);
286  SetHighEnergyLimit(highEnergyLimit[proton]);
287  }
288 
289  if( verboseLevel>0 )
290  {
291  G4cout << "Born ionisation model is initialized " << G4endl
292  << "Energy range: "
293  << LowEnergyLimit() / eV << " eV - "
294  << HighEnergyLimit() / keV << " keV for "
295  << particle->GetParticleName()
296  << G4endl;
297  }
298 
299  // Initialize water density pointer
300  fpMolWaterDensity = G4DNAMolecularMaterial::Instance()->
301  GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER"));
302 
303  //
304  fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
305 
306  if (isInitialised)
307  { return;}
308  fParticleChangeForGamma = GetParticleChangeForGamma();
309  isInitialised = true;
310 }
311 
312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
313 
314 G4double G4DNABornIonisationModel1::CrossSectionPerVolume(const G4Material* material,
315  const G4ParticleDefinition* particleDefinition,
316  G4double ekin,
317  G4double,
318  G4double)
319 {
320  if (verboseLevel > 3)
321  {
322  G4cout << "Calling CrossSectionPerVolume() of G4DNABornIonisationModel1"
323  << G4endl;
324 
325  }
326 
327  if (
328  particleDefinition != G4Proton::ProtonDefinition()
329  &&
330  particleDefinition != G4Electron::ElectronDefinition()
331  )
332 
333  return 0;
334 
335  // Calculate total cross section for model
336 
337  G4double lowLim = 0;
338  G4double highLim = 0;
339  G4double sigma=0;
340 
341  G4double waterDensity = (*fpMolWaterDensity)[material->GetIndex()];
342 
343  if(waterDensity!= 0.0)
344  // if (material == nistwater || material->GetBaseMaterial() == nistwater)
345  {
346  const G4String& particleName = particleDefinition->GetParticleName();
347 
348  std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
349  pos1 = lowEnergyLimit.find(particleName);
350  if (pos1 != lowEnergyLimit.end())
351  {
352  lowLim = pos1->second;
353  }
354 
355  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
356  pos2 = highEnergyLimit.find(particleName);
357  if (pos2 != highEnergyLimit.end())
358  {
359  highLim = pos2->second;
360  }
361 
362  if (ekin >= lowLim && ekin < highLim)
363  {
364  std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
365  pos = tableData.find(particleName);
366 
367  if (pos != tableData.end())
368  {
369  G4DNACrossSectionDataSet* table = pos->second;
370  if (table != 0)
371  {
372  sigma = table->FindValue(ekin);
373 
374  // ICRU49 electronic SP scaling - ZF, SI
375 
376  if (particleDefinition == G4Proton::ProtonDefinition() && ekin < 70*MeV && spScaling)
377  {
378  G4double A = 1.39241700556072800000E-009 ;
379  G4double B = -8.52610412942622630000E-002 ;
380  sigma = sigma * std::exp(A*(ekin/eV)+B);
381  }
382  //
383 
384  }
385  }
386  else
387  {
388  G4Exception("G4DNABornIonisationModel1::CrossSectionPerVolume","em0002",
389  FatalException,"Model not applicable to particle type.");
390  }
391  }
392 
393  if (verboseLevel > 2)
394  {
395  G4cout << "__________________________________" << G4endl;
396  G4cout << "G4DNABornIonisationModel1 - XS INFO START" << G4endl;
397  G4cout << "Kinetic energy(eV)=" << ekin/eV << " particle : " << particleName << G4endl;
398  G4cout << "Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl;
399  G4cout << "Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl;
400  G4cout << "G4DNABornIonisationModel1 - XS INFO END" << G4endl;
401  }
402  } // if (waterMaterial)
403 
404  return sigma*waterDensity;
405 }
406 
407 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
408 
409 void G4DNABornIonisationModel1::SampleSecondaries(std::vector<G4DynamicParticle*>* fvect,
410  const G4MaterialCutsCouple* couple,
411  const G4DynamicParticle* particle,
412  G4double,
413  G4double)
414 {
415 
416  if (verboseLevel > 3)
417  {
418  G4cout << "Calling SampleSecondaries() of G4DNABornIonisationModel1"
419  << G4endl;
420  }
421 
422  G4double lowLim = 0;
423  G4double highLim = 0;
424 
425  G4double k = particle->GetKineticEnergy();
426 
427  const G4String& particleName = particle->GetDefinition()->GetParticleName();
428 
429  std::map< G4String,G4double,std::less<G4String> >::iterator pos1;
430  pos1 = lowEnergyLimit.find(particleName);
431 
432  if (pos1 != lowEnergyLimit.end())
433  {
434  lowLim = pos1->second;
435  }
436 
437  std::map< G4String,G4double,std::less<G4String> >::iterator pos2;
438  pos2 = highEnergyLimit.find(particleName);
439 
440  if (pos2 != highEnergyLimit.end())
441  {
442  highLim = pos2->second;
443  }
444 
445  if (k >= lowLim && k < highLim)
446  {
447  G4ParticleMomentum primaryDirection = particle->GetMomentumDirection();
448  G4double particleMass = particle->GetDefinition()->GetPDGMass();
449  G4double totalEnergy = k + particleMass;
450  G4double pSquare = k * (totalEnergy + particleMass);
451  G4double totalMomentum = std::sqrt(pSquare);
452 
453  G4int ionizationShell = 0;
454 
455  if (!fasterCode) ionizationShell = RandomSelect(k,particleName);
456 
457  // SI: The following protection is necessary to avoid infinite loops :
458  // sigmadiff_ionisation_e_born.dat has non zero partial xs at 18 eV for shell 3 (ionizationShell ==2)
459  // sigmadiff_cumulated_ionisation_e_born.dat has zero cumulated partial xs at 18 eV for shell 3 (ionizationShell ==2)
460  // this is due to the fact that the max allowed transfered energy is (18+10.79)/2=17.025 eV and only transfered energies
461  // strictly above this value have non zero partial xs in sigmadiff_ionisation_e_born.dat (starting at trans = 17.12 eV)
462 
463  if (fasterCode)
464  do
465  {
466  ionizationShell = RandomSelect(k,particleName);
467  }while (k<19*eV && ionizationShell==2 && particle->GetDefinition()==G4Electron::ElectronDefinition());
468 
469  // AM: sample deexcitation
470  // here we assume that H_{2}O electronic levels are the same as Oxygen.
471  // this can be considered true with a rough 10% error in energy on K-shell,
472 
473  G4int secNumberInit = 0;// need to know at a certain point the energy of secondaries
474  G4int secNumberFinal = 0;// So I'll make the diference and then sum the energies
475 
476  G4double bindingEnergy = 0;
477  bindingEnergy = waterStructure.IonisationEnergy(ionizationShell);
478 
479  G4int Z = 8;
480  if(fAtomDeexcitation)
481  {
482  G4AtomicShellEnumerator as = fKShell;
483 
484  if (ionizationShell <5 && ionizationShell >1)
485  {
486  as = G4AtomicShellEnumerator(4-ionizationShell);
487  }
488  else if (ionizationShell <2)
489  {
490  as = G4AtomicShellEnumerator(3);
491  }
492 
493  // FOR DEBUG ONLY
494  // if (ionizationShell == 4) {
495  //
496  // G4cout << "Z: " << Z << " as: " << as
497  // << " ionizationShell: " << ionizationShell << " bindingEnergy: "<< bindingEnergy/eV << G4endl;
498  // G4cout << "Press <Enter> key to continue..." << G4endl;
499  // G4cin.ignore();
500  // }
501 
502  const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
503  secNumberInit = fvect->size();
504  fAtomDeexcitation->GenerateParticles(fvect, shell, Z, 0, 0);
505  secNumberFinal = fvect->size();
506  }
507 
508  G4double secondaryKinetic=-1000*eV;
509 
510  if (fasterCode == false)
511  {
512  secondaryKinetic = RandomizeEjectedElectronEnergy(particle->GetDefinition(),k,ionizationShell);
513  }
514  // SI - 01/04/2014
515  else
516  {
517  secondaryKinetic = RandomizeEjectedElectronEnergyFromCumulatedDcs(particle->GetDefinition(),k,ionizationShell);
518  }
519  //
520 
521  G4ThreeVector deltaDirection =
522  GetAngularDistribution()->SampleDirectionForShell(particle, secondaryKinetic,
523  Z, ionizationShell,
524  couple->GetMaterial());
525 
526  if (particle->GetDefinition() == G4Electron::ElectronDefinition())
527  {
528  G4double deltaTotalMomentum = std::sqrt(secondaryKinetic*(secondaryKinetic + 2.*electron_mass_c2 ));
529 
530  G4double finalPx = totalMomentum*primaryDirection.x() - deltaTotalMomentum*deltaDirection.x();
531  G4double finalPy = totalMomentum*primaryDirection.y() - deltaTotalMomentum*deltaDirection.y();
532  G4double finalPz = totalMomentum*primaryDirection.z() - deltaTotalMomentum*deltaDirection.z();
533  G4double finalMomentum = std::sqrt(finalPx*finalPx + finalPy*finalPy + finalPz*finalPz);
534  finalPx /= finalMomentum;
535  finalPy /= finalMomentum;
536  finalPz /= finalMomentum;
537 
538  G4ThreeVector direction;
539  direction.set(finalPx,finalPy,finalPz);
540 
541  fParticleChangeForGamma->ProposeMomentumDirection(direction.unit());
542  }
543 
544  else fParticleChangeForGamma->ProposeMomentumDirection(primaryDirection);
545 
546  // note that secondaryKinetic is the energy of the delta ray, not of all secondaries.
547  G4double scatteredEnergy = k-bindingEnergy-secondaryKinetic;
548  G4double deexSecEnergy = 0;
549  for (G4int j=secNumberInit; j < secNumberFinal; j++)
550  {
551  deexSecEnergy = deexSecEnergy + (*fvect)[j]->GetKineticEnergy();
552  }
553 
554  if (!statCode)
555  {
556  fParticleChangeForGamma->SetProposedKineticEnergy(scatteredEnergy);
557  fParticleChangeForGamma->ProposeLocalEnergyDeposit(k-scatteredEnergy-secondaryKinetic-deexSecEnergy);
558  }
559  else
560  {
561  fParticleChangeForGamma->SetProposedKineticEnergy(k);
562  fParticleChangeForGamma->ProposeLocalEnergyDeposit(k-scatteredEnergy);
563  }
564 
565  // SI - 01/04/2014
566  if (secondaryKinetic>0)
567  {
568  G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),deltaDirection,secondaryKinetic);
569  fvect->push_back(dp);
570  }
571  //
572 
573  const G4Track * theIncomingTrack = fParticleChangeForGamma->GetCurrentTrack();
574  G4DNAChemistryManager::Instance()->CreateWaterMolecule(eIonizedMolecule,
575  ionizationShell,
576  theIncomingTrack);
577  }
578 }
579 
580 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
581 
582 G4double G4DNABornIonisationModel1::RandomizeEjectedElectronEnergy(G4ParticleDefinition* particleDefinition,
583  G4double k,
584  G4int shell)
585 {
586  // G4cout << "*** SLOW computation for " << " " << particleDefinition->GetParticleName() << G4endl;
587 
588  if (particleDefinition == G4Electron::ElectronDefinition())
589  {
590  G4double maximumEnergyTransfer = 0.;
591  if ((k + waterStructure.IonisationEnergy(shell)) / 2. > k)
592  maximumEnergyTransfer = k;
593  else
594  maximumEnergyTransfer = (k + waterStructure.IonisationEnergy(shell)) / 2.;
595 
596  // SI : original method
597  /*
598  G4double crossSectionMaximum = 0.;
599  for(G4double value=waterStructure.IonisationEnergy(shell); value<=maximumEnergyTransfer; value+=0.1*eV)
600  {
601  G4double differentialCrossSection = DifferentialCrossSection(particleDefinition, k/eV, value/eV, shell);
602  if(differentialCrossSection >= crossSectionMaximum) crossSectionMaximum = differentialCrossSection;
603  }
604  */
605 
606  // SI : alternative method
607  G4double crossSectionMaximum = 0.;
608 
609  G4double minEnergy = waterStructure.IonisationEnergy(shell);
610  G4double maxEnergy = maximumEnergyTransfer;
611  G4int nEnergySteps = 50;
612 
613  G4double value(minEnergy);
614  G4double stpEnergy(std::pow(maxEnergy / value,
615  1. / static_cast<G4double>(nEnergySteps - 1)));
616  G4int step(nEnergySteps);
617  while (step > 0)
618  {
619  step--;
620  G4double differentialCrossSection =
621  DifferentialCrossSection(particleDefinition,
622  k / eV,
623  value / eV,
624  shell);
625  if (differentialCrossSection >= crossSectionMaximum)
626  crossSectionMaximum = differentialCrossSection;
627  value *= stpEnergy;
628  }
629  //
630 
631  G4double secondaryElectronKineticEnergy = 0.;
632  do
633  {
634  secondaryElectronKineticEnergy = G4UniformRand()* (maximumEnergyTransfer-waterStructure.IonisationEnergy(shell));
635  }while(G4UniformRand()*crossSectionMaximum >
636  DifferentialCrossSection(particleDefinition, k/eV,
637  (secondaryElectronKineticEnergy+waterStructure.IonisationEnergy(shell))/eV,shell));
638 
639  return secondaryElectronKineticEnergy;
640 
641  }
642 
643  else if (particleDefinition == G4Proton::ProtonDefinition())
644  {
645  G4double maximumKineticEnergyTransfer = 4.
646  * (electron_mass_c2 / proton_mass_c2) * k;
647 
648  G4double crossSectionMaximum = 0.;
649  for (G4double value = waterStructure.IonisationEnergy(shell);
650  value <= 4. * waterStructure.IonisationEnergy(shell); value += 0.1 * eV)
651  {
652  G4double differentialCrossSection =
653  DifferentialCrossSection(particleDefinition,
654  k / eV,
655  value / eV,
656  shell);
657  if (differentialCrossSection >= crossSectionMaximum)
658  crossSectionMaximum = differentialCrossSection;
659  }
660 
661  G4double secondaryElectronKineticEnergy = 0.;
662  do
663  {
664  secondaryElectronKineticEnergy = G4UniformRand()* maximumKineticEnergyTransfer;
665  }while(G4UniformRand()*crossSectionMaximum >=
666  DifferentialCrossSection(particleDefinition, k/eV,
667  (secondaryElectronKineticEnergy+waterStructure.IonisationEnergy(shell))/eV,shell));
668 
669  return secondaryElectronKineticEnergy;
670  }
671 
672  return 0;
673 }
674 
675 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
676 
677 // The following section is not used anymore but is kept for memory
678 // GetAngularDistribution()->SampleDirectionForShell is used instead
679 
680 /*
681  void G4DNABornIonisationModel1::RandomizeEjectedElectronDirection(G4ParticleDefinition* particleDefinition,
682  G4double k,
683  G4double secKinetic,
684  G4double & cosTheta,
685  G4double & phi )
686  {
687  if (particleDefinition == G4Electron::ElectronDefinition())
688  {
689  phi = twopi * G4UniformRand();
690  if (secKinetic < 50.*eV) cosTheta = (2.*G4UniformRand())-1.;
691  else if (secKinetic <= 200.*eV)
692  {
693  if (G4UniformRand() <= 0.1) cosTheta = (2.*G4UniformRand())-1.;
694  else cosTheta = G4UniformRand()*(std::sqrt(2.)/2);
695  }
696  else
697  {
698  G4double sin2O = (1.-secKinetic/k) / (1.+secKinetic/(2.*electron_mass_c2));
699  cosTheta = std::sqrt(1.-sin2O);
700  }
701  }
702 
703  else if (particleDefinition == G4Proton::ProtonDefinition())
704  {
705  G4double maxSecKinetic = 4.* (electron_mass_c2 / proton_mass_c2) * k;
706  phi = twopi * G4UniformRand();
707 
708  // cosTheta = std::sqrt(secKinetic / maxSecKinetic);
709 
710  // Restriction below 100 eV from Emfietzoglou (2000)
711 
712  if (secKinetic>100*eV) cosTheta = std::sqrt(secKinetic / maxSecKinetic);
713  else cosTheta = (2.*G4UniformRand())-1.;
714 
715  }
716  }
717  */
718 
719 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
720 double G4DNABornIonisationModel1::DifferentialCrossSection(G4ParticleDefinition * particleDefinition,
721  G4double k,
722  G4double energyTransfer,
723  G4int ionizationLevelIndex)
724 {
725  G4double sigma = 0.;
726 
727  if (energyTransfer >= waterStructure.IonisationEnergy(ionizationLevelIndex))
728  {
729  G4double valueT1 = 0;
730  G4double valueT2 = 0;
731  G4double valueE21 = 0;
732  G4double valueE22 = 0;
733  G4double valueE12 = 0;
734  G4double valueE11 = 0;
735 
736  G4double xs11 = 0;
737  G4double xs12 = 0;
738  G4double xs21 = 0;
739  G4double xs22 = 0;
740 
741  if (particleDefinition == G4Electron::ElectronDefinition())
742  {
743  // k should be in eV and energy transfer eV also
744 
745  std::vector<double>::iterator t2 = std::upper_bound(eTdummyVec.begin(),
746  eTdummyVec.end(),
747  k);
748 
749  std::vector<double>::iterator t1 = t2 - 1;
750 
751  // SI : the following condition avoids situations where energyTransfer >last vector element
752  if (energyTransfer <= eVecm[(*t1)].back()
753  && energyTransfer <= eVecm[(*t2)].back())
754  {
755  std::vector<double>::iterator e12 =
756  std::upper_bound(eVecm[(*t1)].begin(),
757  eVecm[(*t1)].end(),
758  energyTransfer);
759  std::vector<double>::iterator e11 = e12 - 1;
760 
761  std::vector<double>::iterator e22 =
762  std::upper_bound(eVecm[(*t2)].begin(),
763  eVecm[(*t2)].end(),
764  energyTransfer);
765  std::vector<double>::iterator e21 = e22 - 1;
766 
767  valueT1 = *t1;
768  valueT2 = *t2;
769  valueE21 = *e21;
770  valueE22 = *e22;
771  valueE12 = *e12;
772  valueE11 = *e11;
773 
774  xs11 = eDiffCrossSectionData[ionizationLevelIndex][valueT1][valueE11];
775  xs12 = eDiffCrossSectionData[ionizationLevelIndex][valueT1][valueE12];
776  xs21 = eDiffCrossSectionData[ionizationLevelIndex][valueT2][valueE21];
777  xs22 = eDiffCrossSectionData[ionizationLevelIndex][valueT2][valueE22];
778 
779  }
780 
781  }
782 
783  if (particleDefinition == G4Proton::ProtonDefinition())
784  {
785  // k should be in eV and energy transfer eV also
786  std::vector<double>::iterator t2 = std::upper_bound(pTdummyVec.begin(),
787  pTdummyVec.end(),
788  k);
789  std::vector<double>::iterator t1 = t2 - 1;
790 
791  std::vector<double>::iterator e12 = std::upper_bound(pVecm[(*t1)].begin(),
792  pVecm[(*t1)].end(),
793  energyTransfer);
794  std::vector<double>::iterator e11 = e12 - 1;
795 
796  std::vector<double>::iterator e22 = std::upper_bound(pVecm[(*t2)].begin(),
797  pVecm[(*t2)].end(),
798  energyTransfer);
799  std::vector<double>::iterator e21 = e22 - 1;
800 
801  valueT1 = *t1;
802  valueT2 = *t2;
803  valueE21 = *e21;
804  valueE22 = *e22;
805  valueE12 = *e12;
806  valueE11 = *e11;
807 
808  xs11 = pDiffCrossSectionData[ionizationLevelIndex][valueT1][valueE11];
809  xs12 = pDiffCrossSectionData[ionizationLevelIndex][valueT1][valueE12];
810  xs21 = pDiffCrossSectionData[ionizationLevelIndex][valueT2][valueE21];
811  xs22 = pDiffCrossSectionData[ionizationLevelIndex][valueT2][valueE22];
812 
813  }
814 
815  G4double xsProduct = xs11 * xs12 * xs21 * xs22;
816  if (xsProduct != 0.)
817  {
818  sigma = QuadInterpolator(valueE11,
819  valueE12,
820  valueE21,
821  valueE22,
822  xs11,
823  xs12,
824  xs21,
825  xs22,
826  valueT1,
827  valueT2,
828  k,
829  energyTransfer);
830  }
831  }
832 
833  return sigma;
834 }
835 
836 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
837 
838 G4double G4DNABornIonisationModel1::Interpolate(G4double e1,
839  G4double e2,
840  G4double e,
841  G4double xs1,
842  G4double xs2)
843 {
844  G4double value = 0.;
845 
846  // Log-log interpolation by default
847 
848  if (e1 != 0 && e2 != 0 && (std::log10(e2) - std::log10(e1)) != 0
849  && !fasterCode)
850  {
851  G4double a = (std::log10(xs2) - std::log10(xs1))
852  / (std::log10(e2) - std::log10(e1));
853  G4double b = std::log10(xs2) - a * std::log10(e2);
854  G4double sigma = a * std::log10(e) + b;
855  value = (std::pow(10., sigma));
856  }
857 
858  // Switch to lin-lin interpolation
859  /*
860  if ((e2-e1)!=0)
861  {
862  G4double d1 = xs1;
863  G4double d2 = xs2;
864  value = (d1 + (d2 - d1)*(e - e1)/ (e2 - e1));
865  }
866  */
867 
868  // Switch to log-lin interpolation for faster code
869  if ((e2 - e1) != 0 && xs1 != 0 && xs2 != 0 && fasterCode)
870  {
871  G4double d1 = std::log10(xs1);
872  G4double d2 = std::log10(xs2);
873  value = std::pow(10., (d1 + (d2 - d1) * (e - e1) / (e2 - e1)));
874  }
875 
876  // Switch to lin-lin interpolation for faster code
877  // in case one of xs1 or xs2 (=cum proba) value is zero
878 
879  if ((e2 - e1) != 0 && (xs1 == 0 || xs2 == 0) && fasterCode)
880  {
881  G4double d1 = xs1;
882  G4double d2 = xs2;
883  value = (d1 + (d2 - d1) * (e - e1) / (e2 - e1));
884  }
885 
886  /*
887  G4cout
888  << e1 << " "
889  << e2 << " "
890  << e << " "
891  << xs1 << " "
892  << xs2 << " "
893  << value
894  << G4endl;
895  */
896 
897  return value;
898 }
899 
900 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
901 
902 G4double G4DNABornIonisationModel1::QuadInterpolator(G4double e11,
903  G4double e12,
904  G4double e21,
905  G4double e22,
906  G4double xs11,
907  G4double xs12,
908  G4double xs21,
909  G4double xs22,
910  G4double t1,
911  G4double t2,
912  G4double t,
913  G4double e)
914 {
915  G4double interpolatedvalue1 = Interpolate(e11, e12, e, xs11, xs12);
916  G4double interpolatedvalue2 = Interpolate(e21, e22, e, xs21, xs22);
917  G4double value = Interpolate(t1,
918  t2,
919  t,
920  interpolatedvalue1,
921  interpolatedvalue2);
922 
923  return value;
924 }
925 
926 G4double G4DNABornIonisationModel1::GetPartialCrossSection(const G4Material* /*material*/,
927  G4int level,
928  const G4ParticleDefinition* particle,
929  G4double kineticEnergy)
930 {
931  std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator pos;
932  pos = tableData.find(particle->GetParticleName());
933 
934  if (pos != tableData.end())
935  {
936  G4DNACrossSectionDataSet* table = pos->second;
937  return table->GetComponent(level)->FindValue(kineticEnergy);
938  }
939 
940  return 0;
941 }
942 
943 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
944 
945 G4int G4DNABornIonisationModel1::RandomSelect(G4double k,
946  const G4String& particle)
947 {
948  G4int level = 0;
949 
950  std::map<G4String, G4DNACrossSectionDataSet*, std::less<G4String> >::iterator pos;
951  pos = tableData.find(particle);
952 
953  if (pos != tableData.end())
954  {
955  G4DNACrossSectionDataSet* table = pos->second;
956 
957  if (table != 0)
958  {
959  G4double* valuesBuffer = new G4double[table->NumberOfComponents()];
960  const size_t n(table->NumberOfComponents());
961  size_t i(n);
962  G4double value = 0.;
963 
964  while (i > 0)
965  {
966  i--;
967  valuesBuffer[i] = table->GetComponent(i)->FindValue(k);
968  value += valuesBuffer[i];
969  }
970 
971  value *= G4UniformRand();
972 
973  i = n;
974 
975  while (i > 0)
976  {
977  i--;
978 
979  if (valuesBuffer[i] > value)
980  {
981  delete[] valuesBuffer;
982  return i;
983  }
984  value -= valuesBuffer[i];
985  }
986 
987  if (valuesBuffer)
988  delete[] valuesBuffer;
989 
990  }
991  } else
992  {
993  G4Exception("G4DNABornIonisationModel1::RandomSelect",
994  "em0002",
995  FatalException,
996  "Model not applicable to particle type.");
997  }
998 
999  return level;
1000 }
1001 
1002 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1003 
1004 G4double G4DNABornIonisationModel1::RandomizeEjectedElectronEnergyFromCumulatedDcs(G4ParticleDefinition* particleDefinition,
1005  G4double k,
1006  G4int shell)
1007 {
1008  //G4cout << "*** FAST computation for " << " " << particleDefinition->GetParticleName() << G4endl;
1009 
1010  G4double secondaryElectronKineticEnergy = 0.;
1011 
1012  G4double random = G4UniformRand();
1013 
1014  secondaryElectronKineticEnergy = TransferedEnergy(particleDefinition,
1015  k / eV,
1016  shell,
1017  random) * eV
1018  - waterStructure.IonisationEnergy(shell);
1019 
1020  //G4cout << RandomTransferedEnergy(particleDefinition, k/eV, shell) << G4endl;
1021  // SI - 01/04/2014
1022  if (secondaryElectronKineticEnergy < 0.)
1023  return 0.;
1024  //
1025 
1026  return secondaryElectronKineticEnergy;
1027 }
1028 
1029 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
1030 
1031 G4double G4DNABornIonisationModel1::TransferedEnergy(G4ParticleDefinition* particleDefinition,
1032  G4double k,
1033  G4int ionizationLevelIndex,
1034  G4double random)
1035 {
1036  G4double nrj = 0.;
1037 
1038  G4double valueK1 = 0;
1039  G4double valueK2 = 0;
1040  G4double valuePROB21 = 0;
1041  G4double valuePROB22 = 0;
1042  G4double valuePROB12 = 0;
1043  G4double valuePROB11 = 0;
1044 
1045  G4double nrjTransf11 = 0;
1046  G4double nrjTransf12 = 0;
1047  G4double nrjTransf21 = 0;
1048  G4double nrjTransf22 = 0;
1049 
1050  if (particleDefinition == G4Electron::ElectronDefinition())
1051  {
1052  // k should be in eV
1053  std::vector<double>::iterator k2 = std::upper_bound(eTdummyVec.begin(),
1054  eTdummyVec.end(),
1055  k);
1056  std::vector<double>::iterator k1 = k2 - 1;
1057 
1058  /*
1059  G4cout << "----> k=" << k
1060  << " " << *k1
1061  << " " << *k2
1062  << " " << random
1063  << " " << ionizationLevelIndex
1064  << " " << eProbaShellMap[ionizationLevelIndex][(*k1)].back()
1065  << " " << eProbaShellMap[ionizationLevelIndex][(*k2)].back()
1066  << G4endl;
1067  */
1068 
1069  // SI : the following condition avoids situations where random >last vector element
1070  if (random <= eProbaShellMap[ionizationLevelIndex][(*k1)].back()
1071  && random <= eProbaShellMap[ionizationLevelIndex][(*k2)].back())
1072  {
1073  std::vector<double>::iterator prob12 =
1074  std::upper_bound(eProbaShellMap[ionizationLevelIndex][(*k1)].begin(),
1075  eProbaShellMap[ionizationLevelIndex][(*k1)].end(),
1076  random);
1077 
1078  std::vector<double>::iterator prob11 = prob12 - 1;
1079 
1080  std::vector<double>::iterator prob22 =
1081  std::upper_bound(eProbaShellMap[ionizationLevelIndex][(*k2)].begin(),
1082  eProbaShellMap[ionizationLevelIndex][(*k2)].end(),
1083  random);
1084 
1085  std::vector<double>::iterator prob21 = prob22 - 1;
1086 
1087  valueK1 = *k1;
1088  valueK2 = *k2;
1089  valuePROB21 = *prob21;
1090  valuePROB22 = *prob22;
1091  valuePROB12 = *prob12;
1092  valuePROB11 = *prob11;
1093 
1094  /*
1095  G4cout << " " << random << " " << valuePROB11 << " "
1096  << valuePROB12 << " " << valuePROB21 << " " << valuePROB22 << G4endl;
1097  */
1098 
1099  nrjTransf11 = eNrjTransfData[ionizationLevelIndex][valueK1][valuePROB11];
1100  nrjTransf12 = eNrjTransfData[ionizationLevelIndex][valueK1][valuePROB12];
1101  nrjTransf21 = eNrjTransfData[ionizationLevelIndex][valueK2][valuePROB21];
1102  nrjTransf22 = eNrjTransfData[ionizationLevelIndex][valueK2][valuePROB22];
1103 
1104  /*
1105  G4cout << " " << ionizationLevelIndex << " "
1106  << random << " " <<valueK1 << " " << valueK2 << G4endl;
1107 
1108  G4cout << " " << random << " " << nrjTransf11 << " "
1109  << nrjTransf12 << " " << nrjTransf21 << " " <<nrjTransf22 << G4endl;
1110  */
1111 
1112  }
1113  // Avoids cases where cum xs is zero for k1 and is not for k2 (with always k1<k2)
1114  if (random > eProbaShellMap[ionizationLevelIndex][(*k1)].back())
1115  {
1116  std::vector<double>::iterator prob22 =
1117  std::upper_bound(eProbaShellMap[ionizationLevelIndex][(*k2)].begin(),
1118  eProbaShellMap[ionizationLevelIndex][(*k2)].end(),
1119  random);
1120 
1121  std::vector<double>::iterator prob21 = prob22 - 1;
1122 
1123  valueK1 = *k1;
1124  valueK2 = *k2;
1125  valuePROB21 = *prob21;
1126  valuePROB22 = *prob22;
1127 
1128  //G4cout << " " << random << " " << valuePROB21 << " " << valuePROB22 << G4endl;
1129 
1130  nrjTransf21 = eNrjTransfData[ionizationLevelIndex][valueK2][valuePROB21];
1131  nrjTransf22 = eNrjTransfData[ionizationLevelIndex][valueK2][valuePROB22];
1132 
1133  G4double interpolatedvalue2 = Interpolate(valuePROB21,
1134  valuePROB22,
1135  random,
1136  nrjTransf21,
1137  nrjTransf22);
1138 
1139  // zeros are explicitely set
1140 
1141  G4double value = Interpolate(valueK1, valueK2, k, 0., interpolatedvalue2);
1142 
1143  /*
1144  G4cout << " " << ionizationLevelIndex << " "
1145  << random << " " <<valueK1 << " " << valueK2 << G4endl;
1146 
1147  G4cout << " " << random << " " << nrjTransf11 << " "
1148  << nrjTransf12 << " " << nrjTransf21 << " " <<nrjTransf22 << G4endl;
1149 
1150  G4cout << "ici" << " " << value << G4endl;
1151  */
1152 
1153  return value;
1154  }
1155  }
1156  //
1157  else if (particleDefinition == G4Proton::ProtonDefinition())
1158  {
1159  // k should be in eV
1160 
1161  std::vector<double>::iterator k2 = std::upper_bound(pTdummyVec.begin(),
1162  pTdummyVec.end(),
1163  k);
1164 
1165  std::vector<double>::iterator k1 = k2 - 1;
1166 
1167  /*
1168  G4cout << "----> k=" << k
1169  << " " << *k1
1170  << " " << *k2
1171  << " " << random
1172  << " " << ionizationLevelIndex
1173  << " " << pProbaShellMap[ionizationLevelIndex][(*k1)].back()
1174  << " " << pProbaShellMap[ionizationLevelIndex][(*k2)].back()
1175  << G4endl;
1176  */
1177 
1178  // SI : the following condition avoids situations where random > last vector element,
1179  // for eg. when the last element is zero
1180  if (random <= pProbaShellMap[ionizationLevelIndex][(*k1)].back()
1181  && random <= pProbaShellMap[ionizationLevelIndex][(*k2)].back())
1182  {
1183  std::vector<double>::iterator prob12 =
1184  std::upper_bound(pProbaShellMap[ionizationLevelIndex][(*k1)].begin(),
1185  pProbaShellMap[ionizationLevelIndex][(*k1)].end(),
1186  random);
1187 
1188  std::vector<double>::iterator prob11 = prob12 - 1;
1189 
1190  std::vector<double>::iterator prob22 =
1191  std::upper_bound(pProbaShellMap[ionizationLevelIndex][(*k2)].begin(),
1192  pProbaShellMap[ionizationLevelIndex][(*k2)].end(),
1193  random);
1194 
1195  std::vector<double>::iterator prob21 = prob22 - 1;
1196 
1197  valueK1 = *k1;
1198  valueK2 = *k2;
1199  valuePROB21 = *prob21;
1200  valuePROB22 = *prob22;
1201  valuePROB12 = *prob12;
1202  valuePROB11 = *prob11;
1203 
1204  /*
1205  G4cout << " " << random << " " << valuePROB11 << " "
1206  << valuePROB12 << " " << valuePROB21 << " " << valuePROB22 << G4endl;
1207  */
1208 
1209  nrjTransf11 = pNrjTransfData[ionizationLevelIndex][valueK1][valuePROB11];
1210  nrjTransf12 = pNrjTransfData[ionizationLevelIndex][valueK1][valuePROB12];
1211  nrjTransf21 = pNrjTransfData[ionizationLevelIndex][valueK2][valuePROB21];
1212  nrjTransf22 = pNrjTransfData[ionizationLevelIndex][valueK2][valuePROB22];
1213 
1214  /*
1215  G4cout << " " << ionizationLevelIndex << " "
1216  << random << " " <<valueK1 << " " << valueK2 << G4endl;
1217 
1218  G4cout << " " << random << " " << nrjTransf11 << " "
1219  << nrjTransf12 << " " << nrjTransf21 << " " <<nrjTransf22 << G4endl;
1220  */
1221  }
1222 
1223  // Avoids cases where cum xs is zero for k1 and is not for k2 (with always k1<k2)
1224 
1225  if (random > pProbaShellMap[ionizationLevelIndex][(*k1)].back())
1226  {
1227  std::vector<double>::iterator prob22 =
1228  std::upper_bound(pProbaShellMap[ionizationLevelIndex][(*k2)].begin(),
1229  pProbaShellMap[ionizationLevelIndex][(*k2)].end(),
1230  random);
1231 
1232  std::vector<double>::iterator prob21 = prob22 - 1;
1233 
1234  valueK1 = *k1;
1235  valueK2 = *k2;
1236  valuePROB21 = *prob21;
1237  valuePROB22 = *prob22;
1238 
1239  //G4cout << " " << random << " " << valuePROB21 << " " << valuePROB22 << G4endl;
1240 
1241  nrjTransf21 = pNrjTransfData[ionizationLevelIndex][valueK2][valuePROB21];
1242  nrjTransf22 = pNrjTransfData[ionizationLevelIndex][valueK2][valuePROB22];
1243 
1244  G4double interpolatedvalue2 = Interpolate(valuePROB21,
1245  valuePROB22,
1246  random,
1247  nrjTransf21,
1248  nrjTransf22);
1249 
1250  // zeros are explicitely set
1251 
1252  G4double value = Interpolate(valueK1, valueK2, k, 0., interpolatedvalue2);
1253 
1254  /*
1255  G4cout << " " << ionizationLevelIndex << " "
1256  << random << " " <<valueK1 << " " << valueK2 << G4endl;
1257 
1258  G4cout << " " << random << " " << nrjTransf11 << " "
1259  << nrjTransf12 << " " << nrjTransf21 << " " <<nrjTransf22 << G4endl;
1260 
1261  G4cout << "ici" << " " << value << G4endl;
1262  */
1263 
1264  return value;
1265  }
1266  }
1267  // End electron and proton cases
1268 
1269  G4double nrjTransfProduct = nrjTransf11 * nrjTransf12 * nrjTransf21
1270  * nrjTransf22;
1271  //G4cout << "nrjTransfProduct=" << nrjTransfProduct << G4endl;
1272 
1273  if (nrjTransfProduct != 0.)
1274  {
1275  nrj = QuadInterpolator(valuePROB11,
1276  valuePROB12,
1277  valuePROB21,
1278  valuePROB22,
1279  nrjTransf11,
1280  nrjTransf12,
1281  nrjTransf21,
1282  nrjTransf22,
1283  valueK1,
1284  valueK2,
1285  k,
1286  random);
1287  }
1288  //G4cout << nrj << endl;
1289 
1290  return nrj;
1291 }
cm
static const double cm
Definition: G4SIunits.hh:118
G4Electron::ElectronDefinition
static G4Electron * ElectronDefinition()
Definition: G4Electron.cc:89
G4DNABornAngle.hh
G4VEmModel::LowEnergyLimit
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:641
MeV
static const double MeV
Definition: G4SIunits.hh:211
G4VEMDataSet::FindValue
virtual G4double FindValue(G4double x, G4int componentId=0) const =0
G4DNABornIonisationModel1::lowEnergyLimit
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
Definition: G4DNABornIonisationModel1.hh:104
d1
static const G4double d1
Definition: G4KleinNishinaCompton.cc:66
G4LossTableManager::Instance
static G4LossTableManager * Instance()
Definition: G4LossTableManager.cc:112
G4DNABornIonisationModel1::fpMolWaterDensity
const std::vector< G4double > * fpMolWaterDensity
Definition: G4DNABornIonisationModel1.hh:99
G4DynamicParticle::GetKineticEnergy
G4double GetKineticEnergy() const
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
G4VEmModel::HighEnergyLimit
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:634
G4DNACrossSectionDataSet::GetComponent
virtual const G4VEMDataSet * GetComponent(G4int componentId) const
Definition: G4DNACrossSectionDataSet.hh:69
G4Material::GetMaterial
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:604
G4DNABornIonisationModel1::RandomizeEjectedElectronEnergyFromCumulatedDcs
G4double RandomizeEjectedElectronEnergyFromCumulatedDcs(G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4int shell)
Definition: G4DNABornIonisationModel1.cc:1004
G4Material::GetIndex
size_t GetIndex() const
Definition: G4Material.hh:262
G4DNABornIonisationModel1::RandomizeEjectedElectronEnergy
G4double RandomizeEjectedElectronEnergy(G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4int shell)
Definition: G4DNABornIonisationModel1.cc:582
G4Proton::ProtonDefinition
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
G4VEmModel::GetAngularDistribution
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:617
e2
static const G4double e2
Definition: G4KleinNishinaCompton.cc:68
G4DNACrossSectionDataSet::LoadData
virtual G4bool LoadData(const G4String &argFileName)
Definition: G4DNACrossSectionDataSet.cc:84
G4LossTableManager.hh
G4DynamicParticle::GetDefinition
G4ParticleDefinition * GetDefinition() const
a
G4double a
Definition: TRTMaterials.hh:39
G4DNABornIonisationModel1.hh
G4DNABornIonisationModel1::fAtomDeexcitation
G4VAtomDeexcitation * fAtomDeexcitation
Definition: G4DNABornIonisationModel1.hh:102
B
double B(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:68
G4DNAMolecularMaterial.hh
G4int
int G4int
Definition: G4Types.hh:78
G4ParticleChangeForGamma::ProposeMomentumDirection
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
Definition: G4ParticleChangeForGamma.hh:147
fKShell
Definition: G4AtomicShellEnumerator.hh:53
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:159
G4VParticleChange::ProposeLocalEnergyDeposit
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4VEmModel::SetHighEnergyLimit
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:725
G4VEmAngularDistribution::SampleDirectionForShell
virtual G4ThreeVector & SampleDirectionForShell(const G4DynamicParticle *dp, G4double finalTotalEnergy, G4int Z, G4int shellID, const G4Material *)
Definition: G4VEmAngularDistribution.cc:60
G4VAtomDeexcitation::GetAtomicShell
virtual const G4AtomicShell * GetAtomicShell(G4int Z, G4AtomicShellEnumerator shell)=0
G4DNABornIonisationModel1::fParticleChangeForGamma
G4ParticleChangeForGamma * fParticleChangeForGamma
Definition: G4DNABornIonisationModel1.hh:90
G4DNABornIonisationModel1::TransferedEnergy
G4double TransferedEnergy(G4ParticleDefinition *aParticleDefinition, G4double incomingParticleEnergy, G4int shell, G4double random)
Definition: G4DNABornIonisationModel1.cc:1031
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:97
G4cout
G4GLOB_DLL std::ostream G4cout
A
double A(double temperature)
Definition: G4DNAElectronHoleRecombination.cc:59
G4DynamicParticle::GetMomentumDirection
const G4ThreeVector & GetMomentumDirection() const
G4DNABornIonisationModel1::QuadInterpolator
G4double QuadInterpolator(G4double e11, G4double e12, G4double e21, G4double e22, G4double x11, G4double x12, G4double x21, G4double x22, G4double t1, G4double t2, G4double t, G4double e)
Definition: G4DNABornIonisationModel1.cc:902
G4DNABornIonisationModel1::highEnergyLimit
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
Definition: G4DNABornIonisationModel1.hh:105
G4DNABornIonisationModel1::Interpolate
G4double Interpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2)
Definition: G4DNABornIonisationModel1.cc:838
G4DNABornIonisationModel1::SampleSecondaries
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
Definition: G4DNABornIonisationModel1.cc:409
n
const G4int n
Definition: G4UrQMD1_3Interface.hh:144
G4DNACrossSectionDataSet::FindValue
virtual G4double FindValue(G4double e, G4int componentId=0) const
Definition: G4DNACrossSectionDataSet.cc:500
G4DNABornIonisationModel1::GetPartialCrossSection
virtual G4double GetPartialCrossSection(const G4Material *, G4int, const G4ParticleDefinition *, G4double)
Definition: G4DNABornIonisationModel1.cc:926
G4DNABornIonisationModel1::waterStructure
G4DNAWaterIonisationStructure waterStructure
Definition: G4DNABornIonisationModel1.hh:124
G4DNACrossSectionDataSet::NumberOfComponents
virtual size_t NumberOfComponents(void) const
Definition: G4DNACrossSectionDataSet.hh:75
e1
static const G4double e1
Definition: G4KleinNishinaCompton.cc:68
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4PhysicalConstants.hh
G4UAtomicDeexcitation.hh
G4DNABornIonisationModel1::Initialise
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &=*(new G4DataVector()))
Definition: G4DNABornIonisationModel1.cc:105
G4DNAChemistryManager.hh
G4DNABornIonisationModel1::G4DNABornIonisationModel1
G4DNABornIonisationModel1(const G4ParticleDefinition *p=0, const G4String &nam="DNABornIonisationModel")
Definition: G4DNABornIonisationModel1.cc:45
G4DNAChemistryManager::Instance
static G4DNAChemistryManager * Instance()
Definition: G4DNAChemistryManager.cc:117
eV
static const double eV
Definition: G4SIunits.hh:212
G4DNAMolecularMaterial::Instance
static G4DNAMolecularMaterial * Instance()
Definition: G4DNAMolecularMaterial.cc:69
eIonizedMolecule
Definition: G4DNAChemistryManager.hh:70
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:161
G4DeltaAngle.hh
G4DNAChemistryManager::CreateWaterMolecule
void CreateWaterMolecule(ElectronicModification, G4int, const G4Track *)
Method used by DNA physics model to create a water molecule.
Definition: G4DNAChemistryManager.cc:572
G4VEmModel::SetAngularDistribution
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:624
G4DNABornIonisationModel1::DifferentialCrossSection
double DifferentialCrossSection(G4ParticleDefinition *aParticleDefinition, G4double k, G4double energyTransfer, G4int shell)
Definition: G4DNABornIonisationModel1.cc:720
FatalException
Definition: G4ExceptionSeverity.hh:60
G4DNABornIonisationModel1::RandomSelect
G4int RandomSelect(G4double energy, const G4String &particle)
Definition: G4DNABornIonisationModel1.cc:945
G4ParticleChangeForGamma::GetCurrentTrack
const G4Track * GetCurrentTrack() const
Definition: G4ParticleChangeForGamma.hh:154
G4Electron::Electron
static G4Electron * Electron()
Definition: G4Electron.cc:94
G4ParticleChangeForGamma::SetProposedKineticEnergy
void SetProposedKineticEnergy(G4double proposedKinEnergy)
Definition: G4ParticleChangeForGamma.hh:129
G4endl
#define G4endl
Definition: G4ios.hh:61
m
static const double m
Definition: G4SIunits.hh:128
G4DNABornIonisationModel1::CrossSectionPerVolume
virtual G4double CrossSectionPerVolume(const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
Definition: G4DNABornIonisationModel1.cc:314
G4LossTableManager::AtomDeexcitation
G4VAtomDeexcitation * AtomDeexcitation()
Definition: G4LossTableManager.hh:530
keV
static const double keV
Definition: G4SIunits.hh:213
G4double
double G4double
Definition: G4Types.hh:76
G4SystemOfUnits.hh
G4VEmModel::SetLowEnergyLimit
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:732
G4VAtomDeexcitation::GenerateParticles
void GenerateParticles(std::vector< G4DynamicParticle * > *secVect, const G4AtomicShell *, G4int Z, G4int coupleIndex)
Definition: G4VAtomDeexcitation.hh:274
G4VEmModel::SetDeexcitationFlag
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:781
d2
static const G4double d2
Definition: G4KleinNishinaCompton.cc:66
G4ParticleMomentum
G4ThreeVector G4ParticleMomentum
Definition: G4ParticleMomentum.hh:50
G4InuclSpecialFunctions::bindingEnergy
G4double bindingEnergy(G4int A, G4int Z)
Definition: bindingEnergy.cc:37
pos
static const G4double pos
Definition: G4ElectroNuclearCrossSection.cc:66
G4AtomicShellEnumerator
G4AtomicShellEnumerator
Definition: G4AtomicShellEnumerator.hh:50
G4VEmModel::GetParticleChangeForGamma
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:134
G4MaterialCutsCouple::GetMaterial
const G4Material * GetMaterial() const
Definition: G4MaterialCutsCouple.hh:150