Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4DNAOneStepThermalizationModel.cc
이 파일의 문서화 페이지로 가기
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: G4DNAOneStepThermalizationModel.cc 110873 2018-06-22 13:11:22Z gcosmo $
27 //
28 // Author: Mathieu Karamitros
29 //
30 // WARNING : This class is released as a prototype.
31 // It might strongly evolve or even disapear in the next releases.
32 //
33 // History:
34 // -----------
35 // 10 Oct 2011 M.Karamitros created
36 //
37 // -------------------------------------------------------------------
38 
39 #include <algorithm>
41 #include "globals.hh"
42 #include "G4Exp.hh"
43 #include "G4RandomDirection.hh"
44 #include "G4Electron.hh"
45 #include "G4EmParameters.hh"
46 
47 //------------------------------------------------------------------------------
48 
49 namespace DNA {
50 namespace Penetration {
51 
52 const double
54 { -4.06217193e-08, 3.06848412e-06, -9.93217814e-05,
55  1.80172797e-03, -2.01135480e-02, 1.42939448e-01,
56  -6.48348714e-01, 1.85227848e+00, -3.36450378e+00,
57  4.37785068e+00, -4.20557339e+00, 3.81679083e+00,
58  -2.34069784e-01 };
59 // fit from Meesungnoen, 2002
60 
61 const double
63 { 0.2, 0.5, 1, 2, 3, 4, 5, 6, 7,
64  // The two last are not in the dataset
65  8, 9}; // eV
66 
67 const double
69 { 17.68*CLHEP::angstrom,
70  22.3*CLHEP::angstrom,
71  28.49*CLHEP::angstrom,
72  45.35*CLHEP::angstrom,
73  70.03*CLHEP::angstrom,
74  98.05*CLHEP::angstrom,
75  120.56*CLHEP::angstrom,
76  132.73*CLHEP::angstrom,
77  142.60*CLHEP::angstrom,
78  // the above value as given in the paper's table does not match
79  // b=27.22 nm nor the mean value. 129.62*CLHEP::angstrom could be
80  // a better fit.
81  //
82  // The two last are made up
83  137.9*CLHEP::angstrom,
84  120.7*CLHEP::angstrom
85 }; // angstrom
86 
87 //----------------------------------------------------------------------------
88 
89 double Meesungnoen2002::GetRmean(double k){
90  G4double k_eV = k/eV;
91 
92  if(k_eV>0.1){ // data until 0.2 eV
93  G4double r_mean = 0;
94  for(int8_t i=12; i!=-1 ; --i){
95  r_mean+=gCoeff[12-i]*std::pow(k_eV,i);
96  }
97  r_mean*=CLHEP::nanometer;
98  return r_mean;
99  }
100  return 0;
101 }
102 
104  G4ThreeVector& displacement)
105 {
106  if(r_mean == 0)
107  {
108  // rare events:
109  // prevent H2O and secondary electron from being placed at the same position
110  displacement = G4RandomDirection() * (1e-3*CLHEP::nanometer);
111  return;
112  }
113 
114  static constexpr double convertRmean3DToSigma1D = 0.62665706865775006;
115  // = sqrt(CLHEP::pi)/pow(2,3./2.)
116 
117  // Use r_mean to build a 3D gaussian
118  const double sigma1D = r_mean * convertRmean3DToSigma1D;
119  displacement = G4ThreeVector(G4RandGauss::shoot(0, sigma1D),
120  G4RandGauss::shoot(0, sigma1D),
121  G4RandGauss::shoot(0, sigma1D));
122 }
123 
125  G4ThreeVector& displacement)
126 {
128 }
129 
130 //----------------------------------------------------------------------------
131 
133  G4ThreeVector& displacement)
134 {
135  GetGaussianPenetrationFromRmean3D(k/eV * 1.8 * nm, // r_mean
136  displacement);
137 }
138 
139 //----------------------------------------------------------------------------
140 
142  G4double k_eV = energy/eV;
143  if(k_eV < 0.2){
144  // rare events:
145  // prevent H2O and secondary electron to be at the spot
146  return 1e-3*CLHEP::nanometer;
147  }
148  else if(k_eV == 9.){
149  return gStdDev_T1990[10];
150  }
151  else if(k_eV > 9.){
152  G4ExceptionDescription description;
153  description << "Terrisol1990 is not tabulated for energies greater than 9eV";
154  G4Exception("Terrisol1990::Get3DStdDeviation",
155  "INVALID_ARGUMENT",
157  description);
158  }
159 
160  size_t lowBin, upBin;
161 
162  if(k_eV >= 1.){
163  lowBin=std::floor(k_eV)+1;
164  upBin=std::min(lowBin+1, size_t(10));
165  }
166  else{
167  auto it=std::lower_bound(&gEnergies_T1990[0],
168  &gEnergies_T1990[2],
169  k_eV);
170  lowBin = it-&gEnergies_T1990[0];
171  upBin = lowBin+1;
172  }
173 
174  double lowE = gEnergies_T1990[lowBin];
175  double upE = gEnergies_T1990[upBin];
176 
177  double lowS = gStdDev_T1990[lowBin];
178  double upS = gStdDev_T1990[upBin];
179 
180  double tanA = (lowS-upS)/(lowE-upE);
181  double sigma3D = lowS + (k_eV-lowE)*tanA;
182  return sigma3D;
183 }
184 
186  double sigma3D=Get3DStdDeviation(energy);
187 
188  static constexpr double s2r=1.595769121605731; // = pow(2,3./2.)/sqrt(CLHEP::pi)
189 
190  double r_mean=sigma3D*s2r;
191  return r_mean;
192 }
193 
195  G4ThreeVector& displacement){
196  double sigma3D = Get3DStdDeviation(energy);
197 
198  static constexpr double factor = 2.20496999539; // = 1./(3. - 8./CLHEP::pi);
199 
200  double sigma1D = std::sqrt(std::pow(sigma3D, 2.)*factor);
201 
202  displacement = G4ThreeVector(G4RandGauss::shoot(0, sigma1D),
203  G4RandGauss::shoot(0, sigma1D),
204  G4RandGauss::shoot(0, sigma1D));
205 }
206 
207 } // Penetration
208 } // DNA
209 
210 //------------------------------------------------------------------------------
211 
213 {
214  G4String modelNamePrefix("DNAOneStepThermalizationModel_");
215 
216  if(penetrationModel == "Terrisol1990")
217  {
218  return new G4TDNAOneStepThermalizationModel<DNA::Penetration::Terrisol1990>(G4Electron::Definition(), modelNamePrefix + penetrationModel);
219  }
220  else if(penetrationModel == "Meesungnoen2002")
221  {
223  }
224  else if(penetrationModel == "Ritchie1994")
225  {
226  return new G4TDNAOneStepThermalizationModel<DNA::Penetration::Ritchie1994>(G4Electron::Definition(), modelNamePrefix + penetrationModel);
227  }
228  else
229  {
230  G4ExceptionDescription description;
231  description << penetrationModel + " is not a valid model name.";
232  G4Exception("G4DNASolvationModelFactory::Create",
233  "INVALID_ARGUMENT",
235  description,
236  "Options are: Terrisol1990, Meesungnoen2002, Ritchie1994.");
237  }
238  return nullptr;
239 }
240 
241 //------------------------------------------------------------------------------
243 {
244  auto dnaSubType = G4EmParameters::Instance()->DNAeSolvationSubType();
245 
246  switch(dnaSubType)
247  {
249  return Create("Ritchie1994");
251  return Create("Terrisol1990");
253  case fDNAUnknownModel:
254  return Create("Meesungnoen2002");
255  default:
256  G4Exception("G4DNASolvationModelFactory::GetMacroDefinedModel",
257  "DnaSubType",
259  "The solvation parameter stored in G4EmParameters is unknown. Supported types are: fRitchie1994eSolvation, fTerrisol1990eSolvation, fMeesungnoen2002eSolvation.");
260  }
261 
262  return nullptr;
263 }
static double Get3DStdDeviation(double energy)
G4DNAModelSubType DNAeSolvationSubType() const
static constexpr double nanometer
Definition: SystemOfUnits.h:81
CLHEP::Hep3Vector G4ThreeVector
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:45
static constexpr double nm
Definition: G4SIunits.hh:112
G4ThreeVector G4RandomDirection()
double G4double
Definition: G4Types.hh:76
static constexpr double angstrom
Definition: SystemOfUnits.h:82
void GetGaussianPenetrationFromRmean3D(G4double r_mean, G4ThreeVector &displacement)
double energy
Definition: plottest35.C:25
static void GetPenetration(G4double energy, G4ThreeVector &displacement)
static constexpr double eV
Definition: G4SIunits.hh:215
static G4VEmModel * Create(const G4String &penetrationModel)
ThreeVector shoot(const G4int Ap, const G4int Af)
static G4Electron * Definition()
Definition: G4Electron.cc:49
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
static void GetPenetration(G4double energy, G4ThreeVector &displacement)
static G4VEmModel * GetMacroDefinedModel()
One step thermalization model can be chosen via macro using /process/dna/e-SolvationSubType Ritchie19...
static G4EmParameters * Instance()
static void GetPenetration(G4double energy, G4ThreeVector &displacement)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments