| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // Author: Alexei Sytov
27
28 /// \file G4ChannelingFastSimInterpolation.cc
29 /// \brief Implementation of the G4ChannelingF
30
31 #include "G4ChannelingFastSimInterpolation.hh"
32 #include "G4SystemOfUnits.hh"
33
34 G4ChannelingFastSimInterpolation::G4Channeling
35
36
37
38
39 {
40 fDx = dx0;
41 fDy = dy0;
42 nPointsx = nPointsx0;
43 nPointsy = nPointsy0;
44 fStepi = fDx/nPointsx;
45 fStepi2= fStepi*fStepi;
46 iModel = iModel0;
47
48 //resize vectors for interpolation coeffic
49 if(iModel==1)
50 {
51 fAI.resize(nPointsx+1);
52 fBI.resize(nPointsx+1);
53 fCI.resize(nPointsx+1);
54 fDI.resize(nPointsx+1);
55 }
56 else if(iModel==2)
57 {
58 fStepj = fDy/nPointsy;
59 fAI3D.resize(nPointsx+1, std::vector<G
60 fBI3D.resize(nPointsx+1, std::vector<G
61 fCI3D.resize(nPointsx+1, std::vector<G
62 fAI3D3.resize(nPointsx+1, std::vector<
63 fBI3D3.resize(nPointsx+1, std::vector<
64 fCI3D3.resize(nPointsx+1, std::vector<
65
66 for(G4int i=0; i<=nPointsx; i++)
67 {
68 fCI3D[i][nPointsy]=0.;
69 }
70 }
71 }
72
73 //....oooOO0OOooo........oooOO0OOooo........oo
74
75 G4double G4ChannelingFastSimInterpolation::Get
76 G4double SplineF=0.;
77 if(iModel==1)
78 {
79 SplineF = Spline1D(xx);
80 }
81 else if(iModel==2)
82 {
83 SplineF = Spline2D(xx, yy);
84 }
85 return SplineF;
86 }
87
88 //....oooOO0OOooo........oooOO0OOooo........oo
89
90 void G4ChannelingFastSimInterpolation::SetCoef
91
92
93
94
95 fAI[i] = AI0;
96 fBI[i] = BI0/cm;
97 fCI[i] = CI0/cm2;
98 fDI[i] = DI0/cm3;
99 }
100
101 //....oooOO0OOooo........oooOO0OOooo........oo
102
103 void G4ChannelingFastSimInterpolation::SetCoef
104
105
106
107
108
109 if (k==0)
110 {
111 fAI3D[i][j] = AI3D0/fStepj/fStepi/6.;
112 fBI3D[i][j] = BI3D0/fStepj/fStepi/6.;
113 fCI3D[i][j] = CI3D0/fStepj/fStepi/6./c
114 }
115 else if (k==1)
116 {
117 fAI3D3[i][j] = AI3D0/fStepj/fStepi/6./
118 fBI3D3[i][j] = BI3D0/fStepj/fStepi/6./
119 fCI3D3[i][j] = CI3D0/fStepj/fStepi/6./
120 }
121 }
122
123 //....oooOO0OOooo........oooOO0OOooo........oo
124
125 G4double G4ChannelingFastSimInterpolation::Spl
126
127 {
128 G4double x1 = xx;
129
130 //if a particle escapes the interpolation
131 if (x1<0.)
132 {
133 x1 += fDx;
134 }
135 else if (x1>=fDx)
136 {
137 x1 -= fDx;
138 }
139
140 // calculation of interpolation function
141 G4int mmx = std::floor(x1/fStepi);
142 x1 -= (mmx+1)*fStepi;
143 G4double Spline3 = fAI[mmx]+x1*(fBI[mmx]+
144 return Spline3;
145 }
146
147 //....oooOO0OOooo........oooOO0OOooo........oo
148
149 G4double G4ChannelingFastSimInterpolation::
150
151 {
152 G4double x1 = xx;
153 G4double y1 = yy;
154
155 //if a particle escapes the interpolation
156 if (x1<0.)
157 {
158 x1 += fDx;
159 }
160 else if (x1>=fDx)
161 {
162 x1 -= fDx;
163 }
164 if (y1<0.)
165 {
166 y1 += fDy;
167 }
168 else if (y1>=fDy)
169 {
170 y1 -= fDy;
171 }
172
173 // calculation of interpolation function
174 G4int mmx = std::floor(x1/fStepi);
175 G4int mmy = std::floor(y1/fStepj);
176
177 G4double tt1 = x1-mmx*fStepi;
178 G4double tt2 = fStepi-tt1;
179 G4double tt13 = tt1*tt1*tt1;
180 G4double tt23 = tt2*tt2*tt2;
181
182 G4double tt1y = y1-mmy*fStepj;
183 G4double tt2y = fStepj-tt1y;
184 G4double tt1y3 = tt1y*tt1y*tt1y;
185 G4double tt2y3 = tt2y*tt2y*tt2y;
186
187 G4double spl3dxx1 = fCI3D3[mmx ][mmy]*tt
188 fAI3D3[mmx ][mmy]*tt
189 G4double spl3dxx2 = fCI3D3[mmx+1][mmy]*tt
190 fAI3D3[mmx+1][mmy]*tt
191 G4double spl3d1 = fCI3D[ mmx ][mmy]*tt
192 fAI3D[ mmx ][mmy]*tt
193 G4double spl3d2 = fCI3D[ mmx+1][mmy]*tt
194 fAI3D[ mmx+1][mmy]*tt
195
196 G4double spl3d = spl3dxx1*tt23 + spl3dxx2
197 (spl3d1*6.-spl3dxx1*fStepi2)*tt2
198
199 return spl3d;
200 }
201
202
203