Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DNA::Penetration::Meesungnoen2002 Struct Reference

#include <G4DNAOneStepThermalizationModel.hh>

Static Public Member Functions

static void GetPenetration (G4double energy, G4ThreeVector &displacement)
 
static double GetRmean (double energy)
 

Static Public Attributes

static const double gCoeff [13]
 

Detailed Description

Definition at line 64 of file G4DNAOneStepThermalizationModel.hh.

Member Function Documentation

void DNA::Penetration::Meesungnoen2002::GetPenetration ( G4double  energy,
G4ThreeVector displacement 
)
static

Definition at line 99 of file G4DNAOneStepThermalizationModel.cc.

100  {
101  displacement=G4ThreeVector(0,0,0);
102  G4double k_eV = k/eV;
103 
104  if(k_eV>0.1){ // data until 0.2 eV
105  G4double r_mean = 0;
106  for(int8_t i=12; i!=-1 ; --i){
107  r_mean+=gCoeff[12-i]*std::pow(k_eV,i);
108  }
109 
110  r_mean*=nanometer;
111 
112  //G4cout << "rmean = " << r_mean << G4endl;
113 
114  static constexpr double r2s=0.62665706865775006; //sqrt(CLHEP::pi)/pow(2,3./2.)
115 
116  // Use r_mean to build a 3D gaussian
117  double sigma3D = r_mean*r2s;
118  double x = G4RandGauss::shoot(0,sigma3D);
119  double y = G4RandGauss::shoot(0,sigma3D);
120  double z = G4RandGauss::shoot(0,sigma3D);
121  displacement=G4ThreeVector(x,y,z);
122  }
123  else{
124  displacement=G4RandomDirection()*(1e-3*CLHEP::nanometer);
125  // rare events:
126  // prevent H2O and secondary electron to be at the spot
127  }
128  }
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4ThreeVector G4RandomDirection()
static constexpr double nanometer
Definition: G4SIunits.hh:101
tuple x
Definition: test.py:50
static constexpr double eV
Definition: G4SIunits.hh:215
static constexpr double nanometer
Definition: SystemOfUnits.h:81
tuple z
Definition: test.py:28
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

double DNA::Penetration::Meesungnoen2002::GetRmean ( double  energy)
static

Definition at line 85 of file G4DNAOneStepThermalizationModel.cc.

85  {
86  G4double k_eV = k/eV;
87 
88  if(k_eV>0.1){ // data until 0.2 eV
89  G4double r_mean = 0;
90  for(int8_t i=12; i!=-1 ; --i){
91  r_mean+=gCoeff[12-i]*std::pow(k_eV,i);
92  }
93  r_mean*=CLHEP::nanometer;
94  return r_mean;
95  }
96  return 0;
97  }
static constexpr double eV
Definition: G4SIunits.hh:215
static constexpr double nanometer
Definition: SystemOfUnits.h:81
double G4double
Definition: G4Types.hh:76

Member Data Documentation

const double DNA::Penetration::Meesungnoen2002::gCoeff
static
Initial value:
=
{ -4.06217193e-08, 3.06848412e-06, -9.93217814e-05,
1.80172797e-03, -2.01135480e-02, 1.42939448e-01,
-6.48348714e-01, 1.85227848e+00, -3.36450378e+00,
4.37785068e+00, -4.20557339e+00, 3.81679083e+00,
-2.34069784e-01 }

Definition at line 70 of file G4DNAOneStepThermalizationModel.hh.


The documentation for this struct was generated from the following files: