Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4GoudsmitSaundersonTable.hh
이 파일의 문서화 페이지로 가기
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: G4GoudsmitSaundersonTable.hh 107824 2017-12-05 15:47:44Z gunter $
27 //
28 // -----------------------------------------------------------------------------
29 //
30 // GEANT4 Class header file
31 //
32 // File name: G4GoudsmitSaundersonTable
33 //
34 // Author: Mihaly Novak / (Omrane Kadri)
35 //
36 // Creation date: 20.02.2009
37 //
38 // Class description:
39 // Class to handle multiple scattering angular distributions precomputed by
40 // using Kawrakow-Bielajew Goudsmit-Saunderson MSC model based on the screened
41 // Rutherford DCS for elastic scattering of electrons/positrons [1,2]. This
42 // class is used by G4GoudsmitSaundersonMscModel to sample the angular
43 // deflection of electrons/positrons after travelling a given path.
44 //
45 // Modifications:
46 // 04.03.2009 V.Ivanchenko cleanup and format according to Geant4 EM style
47 // 18.05.2015 M. Novak This class has been completely replaced (only the original
48 // class name was kept; class description was also inserted):
49 // A new version of Kawrakow-Bielajew Goudsmit-Saunderson MSC model
50 // based on the screened Rutherford DCS for elastic scattering of
51 // electrons/positrons has been introduced[1,2]. The corresponding MSC
52 // angular distributions over a 2D parameter grid have been recomputed
53 // and the CDFs are now stored in a variable transformed (smooth) form
54 // together with the corresponding rational interpolation parameters.
55 // The new version is several times faster, more robust and accurate
56 // compared to the earlier version (G4GoudsmitSaundersonMscModel class
57 // that use these data has been also completely replaced)
58 // 28.04.2017 M. Novak: the GS angular distributions has been recomputed, the
59 // data size has been reduced from 16 MB down to 5 MB by using a new
60 // representation, the class has been modified significantly due to
61 // this new data representation.
62 // 23.08.2017 M. Novak: Added funtionality to handle Mott-correction to the
63 // base GS angular distributions and some other factors (screening
64 // parameter, first and second moments) when Mott-correction is
65 // activated in the GS-MSC model.
66 //
67 // References:
68 // [1] A.F.Bielajew, NIMB, 111 (1996) 195-208
69 // [2] I.Kawrakow, A.F.Bielajew, NIMB 134(1998) 325-336
70 //
71 // -----------------------------------------------------------------------------
72 
73 
74 #ifndef G4GoudsmitSaundersonTable_h
75 #define G4GoudsmitSaundersonTable_h 1
76 
77 #include <vector>
78 
79 #include "G4Types.hh"
80 
81 class G4GSMottCorrection;
83 
85 
86 public:
89 
90  void Initialise(G4double lownergylimit, G4double highenergylimit);
91 
92  // structure to store one GS transformed angular distribution (for a given s/lambda_el,s/lambda_elG1)
93  struct GSMSCAngularDtr {
94  G4int fNumData; // # of data points
95  G4double *fUValues; // array of transformed variables
96  G4double *fParamA; // array of interpolation parameters a
97  G4double *fParamB; // array of interpolation parameters b
98  };
99 
100  void LoadMSCData();
101 
102  G4bool Sampling(G4double lambdaval, G4double qval, G4double scra,
103  G4double &cost, G4double &sint, G4double lekin,
104  G4double beta2, G4int matindx, GSMSCAngularDtr **gsDtr,
105  G4int &mcekini, G4int &mcdelti, G4double &transfPar,
106  G4bool isfirst);
107 
108  G4double SampleCosTheta(G4double lambdaval, G4double qval, G4double scra,
109  G4double lekin, G4double beta2, G4int matindx,
110  GSMSCAngularDtr **gsDtr, G4int &mcekini, G4int &mcdelti,
111  G4double &transfPar, G4bool isfirst);
112 
113  G4double SampleGSSRCosTheta(const GSMSCAngularDtr* gsDrt, G4double transfpar);
114 
115  G4double SingleScattering(G4double lambdaval, G4double scra, G4double lekin,
116  G4double beta2, G4int matindx);
117 
119  G4double &qval, G4double &transfpar);
120 
121  // material dependent MSC parameters (computed at initialisation) regarding
122  // Moliere's screening parameter
123  G4double GetMoliereBc(G4int matindx) { return gMoliereBc[matindx]; }
124 
125  G4double GetMoliereXc2(G4int matindx) { return gMoliereXc2[matindx]; }
126 
127  void GetMottCorrectionFactors(G4double logekin, G4double beta2,
128  G4int matindx, G4double &mcToScr,
129  G4double &mcToQ1, G4double &mcToG2PerG1);
130 
131  // set option to activate/inactivate Mott-correction
133  // set option to activate/inactivate PWA-correction
135 
136  // this method returns with the scattering power correction (to avoid double counting of sub-threshold deflections)
137  // interpolated from tables prepared at initialisation
139 
140  void InitSCPCorrection();
141 
142 private:
143  // initialisation of material dependent Moliere's MSC parameters
144  void InitMoliereMSCParams();
145 
146 
147  private:
148  static G4bool gIsInitialised; // are the precomputed angular distributions already loaded in?
149  static constexpr G4int gLAMBNUM = 64; // # L=s/lambda_el in [fLAMBMIN,fLAMBMAX]
150  static constexpr G4int gQNUM1 = 15; // # Q=s/lambda_el G1 in [fQMIN1,fQMAX1] in the 1-st Q grid
151  static constexpr G4int gQNUM2 = 32; // # Q=s/lambda_el G1 in [fQMIN2,fQMAX2] in the 2-nd Q grid
152  static constexpr G4int gNUMSCR1 = 201; // # of screening parameters in the A(G1) function
153  static constexpr G4int gNUMSCR2 = 51; // # of screening parameters in the A(G1) function
154  static constexpr G4double gLAMBMIN = 1.0; // minimum s/lambda_el
155  static constexpr G4double gLAMBMAX = 100000.0; // maximum s/lambda_el
156  static constexpr G4double gQMIN1 = 0.001; // minimum s/lambda_el G1 in the 1-st Q grid
157  static constexpr G4double gQMAX1 = 0.99; // maximum s/lambda_el G1 in the 1-st Q grid
158  static constexpr G4double gQMIN2 = 0.99; // minimum s/lambda_el G1 in the 2-nd Q grid
159  static constexpr G4double gQMAX2 = 7.99; // maximum s/lambda_el G1 in the 2-nd Q grid
160  //
161  G4bool fIsElectron; // GS-table for e- (for e+ otherwise)
162  G4bool fIsMottCorrection; // flag to indicate if Mott-correction was requested to be used
163  G4bool fIsPWACorrection; // flag to indicate is PWA corrections were requested to be used
164  G4double fLogLambda0; // ln(gLAMBMIN)
165  G4double fLogDeltaLambda; // ln(gLAMBMAX/gLAMBMIN)/(gLAMBNUM-1)
166  G4double fInvLogDeltaLambda; // 1/[ln(gLAMBMAX/gLAMBMIN)/(gLAMBNUM-1)]
167  G4double fInvDeltaQ1; // 1/[(gQMAX1-gQMIN1)/(gQNUM1-1)]
168  G4double fDeltaQ2; // [(gQMAX2-gQMIN2)/(gQNUM2-1)]
169  G4double fInvDeltaQ2; // 1/[(gQMAX2-gQMIN2)/(gQNUM2-1)]
170  //
173  //
174  int fNumSPCEbinPerDec; // scattering power correction energy grid bins per decade
175  struct SCPCorrection {
176  bool fIsUse; //
177  double fPrCut; // sec. e- production cut energy
178  double fLEmin; // log min energy
179  double fILDel; // inverse log delta kinetic energy
180  //std::vector<double> fVEkin; // scattering power correction energies
181  std::vector<double> fVSCPC; // scattering power correction vector
182  };
183  std::vector<SCPCorrection*> fSCPCPerMatCuts;
184 
185 
186  // vector to store all GS transformed angular distributions (cumputed based on the Screened-Rutherford DCS)
187  static std::vector<GSMSCAngularDtr*> gGSMSCAngularDistributions1;
188  static std::vector<GSMSCAngularDtr*> gGSMSCAngularDistributions2;
189 
191 
197  static std::vector<double> gMoliereBc;
198  static std::vector<double> gMoliereXc2;
199  //
200  //
201  G4GSMottCorrection *fMottCorrection;
202 };
203 
204 #endif
std::vector< SCPCorrection * > fSCPCPerMatCuts
void GetMottCorrectionFactors(G4double logekin, G4double beta2, G4int matindx, G4double &mcToScr, G4double &mcToQ1, G4double &mcToG2PerG1)
static constexpr G4double gQMIN1
static constexpr G4double gQMAX1
G4double SingleScattering(G4double lambdaval, G4double scra, G4double lekin, G4double beta2, G4int matindx)
static constexpr G4double gLAMBMIN
G4bool Sampling(G4double lambdaval, G4double qval, G4double scra, G4double &cost, G4double &sint, G4double lekin, G4double beta2, G4int matindx, GSMSCAngularDtr **gsDtr, G4int &mcekini, G4int &mcdelti, G4double &transfPar, G4bool isfirst)
G4GoudsmitSaundersonTable(G4bool iselectron)
void Initialise(G4double lownergylimit, G4double highenergylimit)
G4double GetMoliereBc(G4int matindx)
GSMSCAngularDtr * GetGSAngularDtr(G4double scra, G4double &lambdaval, G4double &qval, G4double &transfpar)
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4double GetMoliereXc2(G4int matindx)
static constexpr G4double gLAMBMAX
G4double SampleGSSRCosTheta(const GSMSCAngularDtr *gsDrt, G4double transfpar)
static std::vector< GSMSCAngularDtr * > gGSMSCAngularDistributions1
G4double ComputeScatteringPowerCorrection(const G4MaterialCutsCouple *matcut, G4double ekin)
static constexpr G4double gQMIN2
int G4int
Definition: G4Types.hh:78
static std::vector< GSMSCAngularDtr * > gGSMSCAngularDistributions2
G4double SampleCosTheta(G4double lambdaval, G4double qval, G4double scra, G4double lekin, G4double beta2, G4int matindx, GSMSCAngularDtr **gsDtr, G4int &mcekini, G4int &mcdelti, G4double &transfPar, G4bool isfirst)
static constexpr G4double gQMAX2