]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | #include "geant321/pilot.h" |
2 | *CMZ : 3.21/02 29/03/94 15.41.24 by S.Giani | |
3 | *-- Author : | |
4 | SUBROUTINE GTGAMA | |
5 | C. | |
6 | C. ****************************************************************** | |
7 | C. * * | |
8 | C. * Photon track. Computes step size and propagates particle * | |
9 | C. * through step. * | |
10 | C. * * | |
11 | C. * ==>Called by : GTRACK * | |
12 | C. * Authors R.Brun, F.Bruyant L.Urban ******** * | |
13 | C. * * | |
14 | C. ****************************************************************** | |
15 | C. | |
16 | #include "geant321/gcbank.inc" | |
17 | #include "geant321/gccuts.inc" | |
18 | #include "geant321/gcjloc.inc" | |
19 | #include "geant321/gconsp.inc" | |
20 | #include "geant321/gcphys.inc" | |
21 | #include "geant321/gcstak.inc" | |
22 | #include "geant321/gctmed.inc" | |
23 | #include "geant321/gcmulo.inc" | |
24 | #include "geant321/gctrak.inc" | |
25 | #if defined(CERNLIB_DEBUG) | |
26 | #include "geant321/gcunit.inc" | |
27 | #endif | |
28 | #if !defined(CERNLIB_SINGLE) | |
29 | PARAMETER (EPSMAC=1.E-6) | |
30 | DOUBLE PRECISION ONE,XCOEF1,XCOEF2,XCOEF3,ZERO | |
31 | #endif | |
32 | #if defined(CERNLIB_SINGLE) | |
33 | PARAMETER (EPSMAC=1.E-11) | |
34 | #endif | |
35 | PARAMETER (ONE=1,ZERO=0) | |
36 | PARAMETER (EPCUT=1.022E-3) | |
37 | C. | |
38 | C. ------------------------------------------------------------------ | |
39 | * | |
40 | * *** Particle below energy threshold ? Short circuit | |
41 | * | |
42 | * | |
43 | IF (GEKIN.LE.CUTGAM) GOTO 998 | |
44 | * | |
45 | * *** Update local pointers if medium has changed | |
46 | * | |
47 | IF(IUPD.EQ.0)THEN | |
48 | IUPD = 1 | |
49 | JPHOT = LQ(JMA-6) | |
50 | JCOMP = LQ(JMA-8) | |
51 | JPAIR = LQ(JMA-10) | |
52 | JPFIS = LQ(JMA-12) | |
53 | JRAYL = LQ(JMA-13) | |
54 | ENDIF | |
55 | * | |
56 | * *** Compute current step size | |
57 | * | |
58 | IPROC = 103 | |
59 | STEP = STEMAX | |
60 | GEKRT1 = 1 .-GEKRAT | |
61 | * | |
62 | * ** Step limitation due to pair production ? | |
63 | * | |
64 | IF (GETOT.GT.EPCUT) THEN | |
65 | IF (IPAIR.GT.0) THEN | |
66 | STEPPA = GEKRT1*Q(JPAIR+IEKBIN) +GEKRAT*Q(JPAIR+IEKBIN+1) | |
67 | SPAIR = STEPPA*ZINTPA | |
68 | IF (SPAIR.LT.STEP) THEN | |
69 | STEP = SPAIR | |
70 | IPROC = 6 | |
71 | ENDIF | |
72 | ENDIF | |
73 | ENDIF | |
74 | * | |
75 | * ** Step limitation due to Compton scattering ? | |
76 | * | |
77 | IF (ICOMP.GT.0) THEN | |
78 | STEPCO = GEKRT1*Q(JCOMP+IEKBIN) +GEKRAT*Q(JCOMP+IEKBIN+1) | |
79 | SCOMP = STEPCO*ZINTCO | |
80 | IF (SCOMP.LT.STEP) THEN | |
81 | STEP = SCOMP | |
82 | IPROC = 7 | |
83 | ENDIF | |
84 | ENDIF | |
85 | * | |
86 | * ** Step limitation due to photo-electric effect ? | |
87 | * | |
88 | IF (GEKIN.LT.0.4) THEN | |
89 | IF (IPHOT.GT.0) THEN | |
90 | STEPPH = GEKRT1*Q(JPHOT+IEKBIN) +GEKRAT*Q(JPHOT+IEKBIN+1) | |
91 | SPHOT = STEPPH*ZINTPH | |
92 | IF (SPHOT.LT.STEP) THEN | |
93 | STEP = SPHOT | |
94 | IPROC = 8 | |
95 | ENDIF | |
96 | ENDIF | |
97 | ENDIF | |
98 | * | |
99 | * ** Step limitation due to photo-fission ? | |
100 | * | |
101 | IF (JPFIS.GT.0) THEN | |
102 | STEPPF = GEKRT1*Q(JPFIS+IEKBIN) +GEKRAT*Q(JPFIS+IEKBIN+1) | |
103 | SPFIS = STEPPF*ZINTPF | |
104 | IF (SPFIS.LT.STEP) THEN | |
105 | STEP = SPFIS | |
106 | IPROC = 23 | |
107 | ENDIF | |
108 | ENDIF | |
109 | * | |
110 | * ** Step limitation due to Rayleigh scattering ? | |
111 | * | |
112 | IF (IRAYL.GT.0) THEN | |
113 | IF (GEKIN.LT.0.01) THEN | |
114 | STEPRA = GEKRT1*Q(JRAYL+IEKBIN) +GEKRAT*Q(JRAYL+IEKBIN+1) | |
115 | SRAYL = STEPRA*ZINTRA | |
116 | IF (SRAYL.LT.STEP) THEN | |
117 | STEP = SRAYL | |
118 | IPROC = 25 | |
119 | ENDIF | |
120 | ENDIF | |
121 | ENDIF | |
122 | * | |
123 | IF (STEP.LT.0.) STEP = 0. | |
124 | * | |
125 | * ** Step limitation due to geometry ? | |
126 | * | |
127 | IF (STEP.GE.SAFETY) THEN | |
128 | CALL GTNEXT | |
129 | IF (IGNEXT.NE.0) THEN | |
130 | STEP = SNEXT + PREC | |
131 | INWVOL= 2 | |
132 | IPROC = 0 | |
133 | NMEC = 1 | |
134 | LMEC(1)=1 | |
135 | ENDIF | |
136 | * | |
137 | * Update SAFETY in stack companions, if any | |
138 | IF (IQ(JSTAK+3).NE.0) THEN | |
139 | DO 10 IST = IQ(JSTAK+3),IQ(JSTAK+1) | |
140 | JST = JSTAK +3 +(IST-1)*NWSTAK | |
141 | Q(JST+11) = SAFETY | |
142 | 10 CONTINUE | |
143 | IQ(JSTAK+3) = 0 | |
144 | ENDIF | |
145 | * | |
146 | ELSE | |
147 | IQ(JSTAK+3) = 0 | |
148 | ENDIF | |
149 | * | |
150 | * *** Linear transport | |
151 | * | |
152 | IF (INWVOL.EQ.2) THEN | |
153 | DO 20 I = 1,3 | |
154 | VECTMP = VECT(I) +STEP*VECT(I+3) | |
155 | IF(VECTMP.EQ.VECT(I)) THEN | |
156 | * | |
157 | * *** Correct for machine precision | |
158 | * | |
159 | IF(VECT(I+3).NE.0.) THEN | |
160 | VECTMP = VECT(I)+ABS(VECT(I))*SIGN(1.,VECT(I+3))* | |
161 | + EPSMAC | |
162 | IF(NMEC.GT.0) THEN | |
163 | IF(LMEC(NMEC).EQ.104) NMEC=NMEC-1 | |
164 | ENDIF | |
165 | NMEC=NMEC+1 | |
166 | LMEC(NMEC)=104 | |
167 | #if defined(CERNLIB_DEBUG) | |
168 | WRITE(CHMAIL, 10000) | |
169 | CALL GMAIL(0,0) | |
170 | WRITE(CHMAIL, 10100) GEKIN, NUMED, STEP, SNEXT | |
171 | CALL GMAIL(0,0) | |
172 | 10000 FORMAT(' Boundary correction in GTGAMA: ', | |
173 | + ' GEKIN NUMED STEP SNEXT') | |
174 | 10100 FORMAT(31X,E10.3,1X,I10,1X,E10.3,1X,E10.3,1X) | |
175 | #endif | |
176 | ENDIF | |
177 | ENDIF | |
178 | VECT(I) = VECTMP | |
179 | 20 CONTINUE | |
180 | ELSE | |
181 | DO 30 I = 1,3 | |
182 | VECT(I) = VECT(I) +STEP*VECT(I+3) | |
183 | 30 CONTINUE | |
184 | ENDIF | |
185 | * | |
186 | SLENG = SLENG +STEP | |
187 | * | |
188 | * *** Update time of flight | |
189 | * | |
190 | TOFG = TOFG +STEP/CLIGHT | |
191 | * | |
192 | * *** Update interaction probabilities | |
193 | * | |
194 | IF (GETOT.GT.EPCUT) THEN | |
195 | IF (IPAIR.GT.0) ZINTPA = ZINTPA -STEP/STEPPA | |
196 | ENDIF | |
197 | IF (ICOMP.GT.0) ZINTCO = ZINTCO -STEP/STEPCO | |
198 | IF (GEKIN.LT.0.4) THEN | |
199 | IF (IPHOT.GT.0) ZINTPH = ZINTPH -STEP/STEPPH | |
200 | ENDIF | |
201 | IF (JPFIS.GT.0) ZINTPF = ZINTPF -STEP/STEPPF | |
202 | IF (IRAYL.GT.0) THEN | |
203 | IF (GEKIN.LT.0.01) ZINTRA = ZINTRA -STEP/STEPRA | |
204 | ENDIF | |
205 | * | |
206 | IF (IPROC.EQ.0) GO TO 999 | |
207 | NMEC = 1 | |
208 | LMEC(1) = IPROC | |
209 | * | |
210 | * ** Pair production ? | |
211 | * | |
212 | IF (IPROC.EQ.6) THEN | |
213 | CALL GPAIRG | |
214 | * | |
215 | * ** Compton scattering ? | |
216 | * | |
217 | ELSE IF (IPROC.EQ.7) THEN | |
218 | CALL GCOMP | |
219 | * | |
220 | * ** Photo-electric effect ? | |
221 | * | |
222 | ELSE IF (IPROC.EQ.8) THEN | |
223 | * Calculate range of the photoelectron ( with kin. energy Ephot) | |
224 | * | |
225 | IF(GEKIN.LE.0.001) THEN | |
226 | JCOEF = LQ(JMA-17) | |
227 | IF(GEKRAT.LT.0.7) THEN | |
228 | I1 = MAX(IEKBIN-1,1) | |
229 | ELSE | |
230 | I1 = MIN(IEKBIN,NEKBIN-1) | |
231 | ENDIF | |
232 | I1 = 3*(I1-1)+1 | |
233 | XCOEF1 = Q(JCOEF+I1) | |
234 | XCOEF2 = Q(JCOEF+I1+1) | |
235 | XCOEF3 = Q(JCOEF+I1+2) | |
236 | IF(XCOEF1.NE.0.) THEN | |
237 | STOPMX = -XCOEF2+SIGN(ONE,XCOEF1)*SQRT(XCOEF2**2 - (XCOEF3- | |
238 | + GEKIN/XCOEF1)) | |
239 | ELSE | |
240 | STOPMX = - (XCOEF3-GEKIN)/XCOEF2 | |
241 | ENDIF | |
242 | * | |
243 | * DO NOT call GPHOT if this (overestimated) range is smaller | |
244 | * than SAFETY | |
245 | * | |
246 | IF (STOPMX.LE.SAFETY) GOTO 998 | |
247 | ENDIF | |
248 | ||
249 | CALL GPHOT | |
250 | * | |
251 | * ** Rayleigh effect ? | |
252 | * | |
253 | ELSE IF (IPROC.EQ.25) THEN | |
254 | CALL GRAYL | |
255 | * | |
256 | * ** Photo-fission ? | |
257 | * | |
258 | ELSE IF (IPROC.EQ.23) THEN | |
259 | CALL GPFIS | |
260 | * | |
261 | ENDIF | |
262 | * | |
263 | GOTO 999 | |
264 | 998 DESTEP = GEKIN | |
265 | GEKIN = 0. | |
266 | GETOT = 0. | |
267 | VECT(7)= 0. | |
268 | ISTOP = 2 | |
269 | NMEC = 1 | |
270 | LMEC(1)= 30 | |
271 | 999 END |