Geant4-11
G4DNARuddIonisationExtendedModel.hh
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//
27
28#ifndef G4DNARuddIonisationExtendedModel_h
29#define G4DNARuddIonisationExtendedModel_h 1
30
31#include "G4VEmModel.hh"
34
37#include "G4Electron.hh"
38#include "G4Proton.hh"
40
43#include "G4NistManager.hh"
44
46{
47
48public:
49
51 const G4String& nam = "DNARuddIonisationExtendedModel");
52
54
55 virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&);
56
58 const G4ParticleDefinition* p,
59 G4double ekin,
60 G4double emin,
62
63 virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
65 const G4DynamicParticle*,
66 G4double tmin,
67 G4double maxEnergy);
68
69 inline void SelectStationary(G4bool input);
70
71protected:
72
74
75private:
76
78
79 // Water density table
80 const std::vector<G4double>* fpWaterDensity;
81
82 //deexcitation manager to produce fluo photns and e-
84
85 std::map<G4String,G4double,std::less<G4String> > lowEnergyLimit;
86 std::map<G4String,G4double,std::less<G4String> > highEnergyLimit;
87
88 // ZF: 26-10-2010
90
93
94 // Cross section
95
96 typedef std::map<G4String,G4String,std::less<G4String> > MapFile;
98
99 typedef std::map<G4String,G4DNACrossSectionDataSet*,std::less<G4String> > MapData;
101
102 // Final state
103
105
107 G4double incomingParticleEnergy,
108 G4int shell);
109
110 // SI: Not necessary anymore since we now use interface to angle generator
111 /*
112 void RandomizeEjectedElectronDirection(G4ParticleDefinition* particleDefinition,
113 G4double incomingParticleEnergy,
114 G4double outgoingParticleEnergy,
115 G4double & cosTheta,
116 G4double & phi, G4int shell);
117 */
118
120 G4double k,
121 G4double proposed_ws,
122 G4int ionizationLevelIndex);
123
125 G4double k,
126 G4int ionizationLevelIndex);
127
128 G4double CorrectionFactor(G4ParticleDefinition* particleDefinition, G4double k, G4int shell);
129
131 G4double energyTransferred,
132 G4double slaterEffectiveChg,
133 G4double shellNumber);
134
136 G4double energyTransferred,
137 G4double slaterEffectiveChg,
138 G4double shellNumber);
139
140
142 G4double energyTransferred,
143 G4double slaterEffectiveChg,
144 G4double shellNumber);
145
147 G4double energyTransferred,
148 G4double slaterEffectiveChg,
149 G4double shellNumber) ;
150
153
154 // Partial cross section
155
157
158 G4double Sum(G4double energy, const G4String& particle);
159
160 G4int RandomSelect(G4double energy,const G4String& particle );
161
162 //
163
166
167};
168
169//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
170
172{
173 statCode = input;
174}
175
176//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
177
178#endif
static const G4double emax
double G4double
Definition: G4Types.hh:83
bool G4bool
Definition: G4Types.hh:86
int G4int
Definition: G4Types.hh:85
G4double CorrectionFactor(G4ParticleDefinition *particleDefinition, G4double k, G4int shell)
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
std::map< G4double, G4double > killBelowEnergyForA
G4double ProposedSampledEnergy(G4ParticleDefinition *particle, G4double k, G4int ionizationLevelIndex)
std::map< G4double, G4double > lowEnergyLimitOfModelForA
G4double R(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
G4DNARuddIonisationExtendedModel(const G4DNARuddIonisationExtendedModel &)
G4double PartialCrossSection(const G4Track &track)
G4double RandomizeEjectedElectronEnergy(G4ParticleDefinition *particleDefinition, G4double incomingParticleEnergy, G4int shell)
G4double RejectionFunction(G4ParticleDefinition *particle, G4double k, G4double proposed_ws, G4int ionizationLevelIndex)
std::map< G4String, G4double, std::less< G4String > > lowEnergyLimit
G4double S_2s(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
std::map< G4String, G4double, std::less< G4String > > highEnergyLimit
virtual G4double CrossSectionPerVolume(const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
std::map< G4double, G4double > lowEnergyLimitForA
G4double S_1s(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
std::map< G4String, G4DNACrossSectionDataSet *, std::less< G4String > > MapData
const std::vector< G4double > * fpWaterDensity
G4double S_2p(G4double t, G4double energyTransferred, G4double slaterEffectiveChg, G4double shellNumber)
std::map< G4String, G4String, std::less< G4String > > MapFile
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)
G4DNARuddIonisationExtendedModel & operator=(const G4DNARuddIonisationExtendedModel &right)
G4DNARuddIonisationExtendedModel(const G4ParticleDefinition *p=0, const G4String &nam="DNARuddIonisationExtendedModel")
G4int RandomSelect(G4double energy, const G4String &particle)
G4double Sum(G4double energy, const G4String &particle)
G4double energy(const ThreeVector &p, const G4double m)
string material
Definition: eplot.py:19