1 | SUBROUTINE PHOTON(IRCODE)
|
---|
2 | C VERSION 4.00 -- 26 JAN 1986/1900
|
---|
3 | C******************************************************************
|
---|
4 | DOUBLE PRECISION PEIG
|
---|
5 | COMMON/BOUNDS/ECUT(6),PCUT(6),VACDST
|
---|
6 | *KEEP,BUFFS.
|
---|
7 | COMMON /BUFFS/ RUNH,RUNE,EVTH,EVTE,DATAB,LH
|
---|
8 | INTEGER MAXBUF,MAXLEN
|
---|
9 | PARAMETER (MAXBUF=39*7)
|
---|
10 | PARAMETER (MAXLEN=12)
|
---|
11 | REAL RUNH(MAXBUF),EVTH(MAXBUF),EVTE(MAXBUF),
|
---|
12 | * RUNE(MAXBUF),DATAB(MAXBUF)
|
---|
13 | INTEGER LH
|
---|
14 | CHARACTER*4 CRUNH,CRUNE,CEVTH,CEVTE
|
---|
15 | EQUIVALENCE (RUNH(1),CRUNH), (RUNE(1),CRUNE)
|
---|
16 | EQUIVALENCE (EVTH(1),CEVTH), (EVTE(1),CEVTE)
|
---|
17 | *KEEP,CEREN1.
|
---|
18 | COMMON /CEREN1/ CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD,
|
---|
19 | * CERSIZ,LCERFI
|
---|
20 | DOUBLE PRECISION CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD
|
---|
21 | REAL CERSIZ
|
---|
22 | LOGICAL LCERFI
|
---|
23 | *KEEP,EPCONT.
|
---|
24 | COMMON/EPCONT/ EDEP,RATIO,TSTEP,TUSTEP,USTEP,TVSTEP,VSTEP,IDISC,
|
---|
25 | * IROLD,IRNEW,RHOFAC, EOLD,ENEW,EKE,ELKE,BETA2,GLE,
|
---|
26 | * TSCAT,IAUSFL
|
---|
27 | DOUBLE PRECISION EDEP,RATIO
|
---|
28 | REAL TSTEP,TUSTEP,USTEP,TVSTEP,VSTEP,RHOFAC,EOLD,ENEW,
|
---|
29 | * EKE,ELKE,BETA2,GLE,TSCAT
|
---|
30 | INTEGER IDISC,IROLD,IRNEW,IAUSFL(29)
|
---|
31 | *KEND.
|
---|
32 | COMMON/GEOM/ZALTIT,BOUND(6),NEWOBS,OBSLVL(10)
|
---|
33 | *KEEP,LONGI.
|
---|
34 | COMMON /LONGI/ APLONG,HLONG,PLONG,SPLONG,THSTEP,THSTPI,
|
---|
35 | * NSTEP,LLONGI,FLGFIT
|
---|
36 | DOUBLE PRECISION APLONG(0:1040,9),HLONG(0:1024),PLONG(0:1040,9),
|
---|
37 | * SPLONG(0:1040,9),THSTEP,THSTPI
|
---|
38 | INTEGER NSTEP
|
---|
39 | LOGICAL LLONGI,FLGFIT
|
---|
40 | *KEND.
|
---|
41 | COMMON /MEDIA/ NMED, RLC,RLDU,RLDUI,RHO,MSGE,MGE,MSEKE,MEKE,MLEKE,
|
---|
42 | *MCMFP,MRANGE,IRAYLM,HBARO(6),HBAROI(6)
|
---|
43 | CHARACTER MEDIA*24
|
---|
44 | COMMON/MEDIAC/MEDIA
|
---|
45 | COMMON/MISC/KMPI,KMPO,DUNIT,NOSCAT,MED(6),RHOR(6),IRAYLR(6)
|
---|
46 | DOUBLE PRECISION PRRMMU
|
---|
47 | COMMON/MUON/PRRMMU,RMMU,RMMUT2
|
---|
48 | *KEEP,OBSPAR.
|
---|
49 | COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP,
|
---|
50 | * THETPR,PHIPR,NOBSLV
|
---|
51 | DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10),
|
---|
52 | * THETAP,THETPR(2),PHIP,PHIPR(2)
|
---|
53 | INTEGER NOBSLV
|
---|
54 | *KEEP,PARPAR.
|
---|
55 | COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
|
---|
56 | * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
|
---|
57 | DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
|
---|
58 | * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
|
---|
59 | INTEGER ITYPE,LEVL
|
---|
60 | *KEND.
|
---|
61 | COMMON/PHOTIN/EBINDA,GE0,GE1, MPGEM(1),GMFP0(500),GMFP1(500),GBR10
|
---|
62 | *(500),GBR11(500),GBR20(500),GBR21(500),GBR30(500),GBR31(500),GBR40
|
---|
63 | *(500),GBR41(500),NGR,RCO0,RCO1, RSCT0(100),RSCT1(100), COHE0(500),
|
---|
64 | *COHE1(500)
|
---|
65 | DOUBLE PRECISION PI0MSQ
|
---|
66 | COMMON/PION/PI0MSQ,PITHR,PICMAS,PI0MAS,AMASK0,AMASKC,AMASPR,AMASNT
|
---|
67 | *
|
---|
68 | *KEEP,RANDPA.
|
---|
69 | COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
|
---|
70 | DOUBLE PRECISION FAC,U1,U2
|
---|
71 | REAL RD(3000)
|
---|
72 | INTEGER ISEED(103,10),NSEQ
|
---|
73 | LOGICAL KNOR
|
---|
74 | *KEEP,REJECT.
|
---|
75 | COMMON /REJECT/ AVNREJ,
|
---|
76 | * ALTMIN,ANEXP,THICKA,THICKD,CUTLN,EONCUT,
|
---|
77 | * FNPRIM
|
---|
78 | DOUBLE PRECISION AVNREJ(10)
|
---|
79 | REAL ALTMIN(10),ANEXP(10),THICKA(10),THICKD(10),
|
---|
80 | * CUTLN,EONCUT
|
---|
81 | LOGICAL FNPRIM
|
---|
82 | *KEEP,RUNPAR.
|
---|
83 | COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
|
---|
84 | * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
|
---|
85 | * MONIOU,MDEBUG,NUCNUC,
|
---|
86 | * CETAPE,
|
---|
87 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
---|
88 | * N1STTR,MDBASE,
|
---|
89 | * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
---|
90 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
---|
91 | * ,GHEISH,GHESIG
|
---|
92 | COMMON /RUNPAC/ DSN,HOST,USER
|
---|
93 | DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
|
---|
94 | REAL STEPFC
|
---|
95 | INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
|
---|
96 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
---|
97 | * N1STTR,MDBASE
|
---|
98 | INTEGER CETAPE
|
---|
99 | CHARACTER*79 DSN
|
---|
100 | CHARACTER*20 HOST,USER
|
---|
101 |
|
---|
102 | LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
---|
103 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
---|
104 | * ,GHEISH,GHESIG
|
---|
105 | *KEEP,STACKE.
|
---|
106 | COMMON/STACKE/ E,TIME,X,Y,Z,U,V,W,DNEAR,IQ,IGEN,IR,IOBS,LPCTE,NP
|
---|
107 | DOUBLE PRECISION E(60),TIME(60)
|
---|
108 | REAL X(60),Y(60),Z(60),U(60),V(60),W(60),DNEAR(60)
|
---|
109 | INTEGER IQ(60),IGEN(60),IR(60),IOBS(60),LPCTE(60),NP
|
---|
110 | *KEND.
|
---|
111 | COMMON/THRESH/RMT2,RMSQ,ESCD2,AP,API,AE,UP,UE,TE,THMOLL
|
---|
112 | COMMON/UPHIOT/THETA,SINTHE,COSTHE,SINPHI, COSPHI,PI,TWOPI,PI5D2
|
---|
113 | DOUBLE PRECISION PZERO,PRM,PRMT2,RMI,VC
|
---|
114 | COMMON/USEFUL/PZERO,PRM,PRMT2,RMI,VC,RM,MEDIUM,MEDOLD,IBLOBE,ICALL
|
---|
115 | COMMON/ACLOCK/NCLOCK,JCLOCK
|
---|
116 | DOUBLE PRECISION THICK
|
---|
117 |
|
---|
118 |
|
---|
119 |
|
---|
120 | C_____IF (NCLOCK.GT.JCLOCK) THEN
|
---|
121 | C______WRITE(MDEBUG,* )' PHOTON:NP=',NP,' IR=',IR(NP),' IOBS=',IOBS(NP)
|
---|
122 | C______CALL AUSGB2
|
---|
123 | C_____END IF
|
---|
124 | NEWOBS=IOBS(NP)
|
---|
125 | IRCODE=1
|
---|
126 | PEIG=E(NP)
|
---|
127 | EIG=PEIG
|
---|
128 | IRL=IR(NP)
|
---|
129 | MEDIUM=MED(IRL)
|
---|
130 | IF((EIG.LE.PCUT(IRL)))GO TO 970
|
---|
131 | 980 CONTINUE
|
---|
132 | 981 CONTINUE
|
---|
133 | GLE=LOG(EIG)
|
---|
134 | CALL RMMAR(RNNO35,1,2)
|
---|
135 | IF ((RNNO35.EQ.0.0)) THEN
|
---|
136 | RNNO35=1.E-30
|
---|
137 | END IF
|
---|
138 | DPMFP=-ALOG(RNNO35)
|
---|
139 | IROLD=IR(NP)
|
---|
140 | 1030 CONTINUE
|
---|
141 | 1031 CONTINUE
|
---|
142 | IF (MEDIUM.NE.0) THEN
|
---|
143 | LGLE=GE1*GLE+GE0
|
---|
144 | GMFPR0=GMFP1(LGLE)*GLE+GMFP0(LGLE)
|
---|
145 | END IF
|
---|
146 | 1040 CONTINUE
|
---|
147 | 1041 CONTINUE
|
---|
148 | IF (MEDIUM.EQ.0) THEN
|
---|
149 | TSTEP=VACDST
|
---|
150 | ELSE
|
---|
151 | RHOFAC=RHOR(IRL)/RHO
|
---|
152 | RHOFI=1./RHOFAC
|
---|
153 | GMFP=GMFPR0*RHOFI
|
---|
154 | IF ((IRAYLR(IRL).EQ.1)) THEN
|
---|
155 | COHFAC=COHE1(LGLE)*GLE+COHE0(LGLE)
|
---|
156 | GMFP=GMFP*COHFAC
|
---|
157 | END IF
|
---|
158 | TSTEP=GMFP*DPMFP
|
---|
159 | ALTEXP=EXP(-Z(NP)*HBAROI(IRL))
|
---|
160 | TSTEP=TSTEP*ALTEXP
|
---|
161 | DISC=W(NP)*TSTEP*HBAROI(IRL)
|
---|
162 | IF (ABS(DISC).LT.0.065) THEN
|
---|
163 | TSTEP=TSTEP*(1.-0.5*DISC*(1.-0.6666667*DISC* (1.-0.75*DISC *
|
---|
164 | * (1.-0.8*DISC))))
|
---|
165 | ELSE IF(DISC.LE.-1.) THEN
|
---|
166 | TSTEP=VACDST
|
---|
167 | ELSE
|
---|
168 | TSTEP=TSTEP/DISC*LOG(DISC+1.)
|
---|
169 | END IF
|
---|
170 | END IF
|
---|
171 | IRNEW=IR(NP)
|
---|
172 | IDISC=0
|
---|
173 | USTEP=TSTEP
|
---|
174 | TUSTEP=USTEP
|
---|
175 | IF((USTEP.GT.DNEAR(NP)))CALL HOWFAR
|
---|
176 | IF((IDISC.GT.0))GO TO 1000
|
---|
177 | VSTEP=USTEP
|
---|
178 | TVSTEP=VSTEP
|
---|
179 | EDEP=PZERO
|
---|
180 | USTEPU=USTEP
|
---|
181 | DISC=W(NP)*USTEPU*HBAROI(IRL)
|
---|
182 | USTEPU=USTEPU/ALTEXP
|
---|
183 | IF (ABS(DISC).LT.0.16) THEN
|
---|
184 | USTEPU=USTEPU*(1.+.5*DISC*(1.+.33333333*DISC* (1.+0.25*DISC* (
|
---|
185 | * 1.+0.2*DISC))))
|
---|
186 | ELSE
|
---|
187 | USTEPU=USTEPU/DISC*(EXP(DISC)-1.)
|
---|
188 | END IF
|
---|
189 | X(NP)=X(NP)+U(NP)*USTEP
|
---|
190 | Y(NP)=Y(NP)+V(NP)*USTEP
|
---|
191 | ZOLD =Z(NP)
|
---|
192 | Z(NP)=Z(NP)+W(NP)*USTEP
|
---|
193 | TIME(NP)=TIME(NP)+TVSTEP*VC
|
---|
194 |
|
---|
195 |
|
---|
196 | C ADD PHOTONS TO THE LONGITUDINAL DEVELOPMENT
|
---|
197 | IF ( LLONGI ) THEN
|
---|
198 | C FIND FIRST THE EQUIVALENT LEVELS
|
---|
199 | C IF STARTING POINT BELOW LOWEST LEVEL THEN DON'T CHECK
|
---|
200 | IF ( HLONG(NSTEP) .LE. -ZOLD ) THEN
|
---|
201 | LPCT1 = LPCTE(NP)
|
---|
202 | C Z NEW IS PROBABLY ONLY LITTLE BELOW Z OLD, THEREFORE INCREMENTAL SEARCH
|
---|
203 | DO 6002 I1 = LPCT1,NSTEP
|
---|
204 | IF ( HLONG(I1) .LT. -Z(NP) ) GOTO 6003
|
---|
205 | 6002 CONTINUE
|
---|
206 | I1 = NSTEP + 1
|
---|
207 | 6003 CONTINUE
|
---|
208 | LPCT2 = I1 - 1
|
---|
209 | DO 485 I=LPCT1,LPCT2
|
---|
210 | PLONG(I,1) = PLONG(I,1) + 1.D0
|
---|
211 | 485 CONTINUE
|
---|
212 | LPCTE(NP) = LPCT2 + 1
|
---|
213 | ENDIF
|
---|
214 | ENDIF
|
---|
215 | DNEAR(NP)=DNEAR(NP)-USTEP
|
---|
216 | IF (MEDIUM.NE.0) THEN
|
---|
217 | DPMFP=MAX(0.,DPMFP-USTEPU/GMFP)
|
---|
218 | END IF
|
---|
219 | IROLD=IR(NP)
|
---|
220 | MEDOLD=MEDIUM
|
---|
221 | IF (IRNEW.NE.IROLD) THEN
|
---|
222 | IR(NP)=IRNEW
|
---|
223 | IRL=IRNEW
|
---|
224 | MEDIUM=MED(IRL)
|
---|
225 | IF((EIG.LE.PCUT(IRL)))GO TO 970
|
---|
226 | END IF
|
---|
227 | IF (NEWOBS.GT.IOBS(NP)) THEN
|
---|
228 | CALL AUSGAB
|
---|
229 | IOBS(NP)=NEWOBS
|
---|
230 | END IF
|
---|
231 | IF((IDISC.LT.0))GO TO 1000
|
---|
232 | IF((MEDIUM.NE.MEDOLD))GO TO 1042
|
---|
233 | IF((MEDIUM.NE.0.AND.DPMFP.LE.1.E-6))GO TO 1032
|
---|
234 | GO TO 1041
|
---|
235 | 1042 CONTINUE
|
---|
236 | GO TO 1031
|
---|
237 | 1032 CONTINUE
|
---|
238 | IF ((IRAYLR(IRL).EQ.1)) THEN
|
---|
239 | CALL RMMAR(RNNO37,1,2)
|
---|
240 | IF ((RNNO37.LE.(1.0-COHFAC))) THEN
|
---|
241 | 1050 CONTINUE
|
---|
242 | 1051 CONTINUE
|
---|
243 | CALL RMMAR(XXX,1,2)
|
---|
244 | LXXX=RCO1*XXX+RCO0
|
---|
245 | X2=RSCT1(LXXX)*XXX+RSCT0(LXXX)
|
---|
246 | Q2=X2*RMSQ*.23547885E-02
|
---|
247 | COSTHE=1.-Q2/(2.*E(NP)*E(NP))
|
---|
248 | IF((ABS(COSTHE).GT.1.0))GO TO 1050
|
---|
249 | CSQTHE=COSTHE*COSTHE
|
---|
250 | REJF=(1.0+CSQTHE)*.5
|
---|
251 | CALL RMMAR(RNNORJ,1,2)
|
---|
252 | IF((RNNORJ.LE.REJF))GO TO1052
|
---|
253 | GO TO 1051
|
---|
254 | 1052 CONTINUE
|
---|
255 | SINTHE=SQRT(AMAX1(0.,1.0-CSQTHE))
|
---|
256 | CALL UPHI(2,1)
|
---|
257 | GOTO 980
|
---|
258 | END IF
|
---|
259 | END IF
|
---|
260 | IF ( .NOT. FNPRIM ) THEN
|
---|
261 | X(1)=0.
|
---|
262 | Y(1)=0.
|
---|
263 | EVTH(5)=X(1)
|
---|
264 | EVTH(6)=-Y(1)
|
---|
265 | IF (FIX1I) THEN
|
---|
266 | Z(1)=-FIXHEI
|
---|
267 | NP=1
|
---|
268 | LPCTE(1)=MIN(NSTEP,INT(THICK(FIXHEI)*THSTPI)+1)
|
---|
269 | SITHET=SQRT(1.D0-SECPAR(3)**2)
|
---|
270 | U(1)=SITHET*COS(-SECPAR(4))
|
---|
271 | V(1)=SITHET*SIN(-SECPAR(4))
|
---|
272 | W(1)=SECPAR(3)
|
---|
273 | RADINV=1.5-0.5*(U(1)**2+V(1)**2+W(1)**2)
|
---|
274 | U(1)=U(1)*RADINV
|
---|
275 | V(1)=V(1)*RADINV
|
---|
276 | W(1)=W(1)*RADINV
|
---|
277 | END IF
|
---|
278 | EVTH(7)=-Z(1)
|
---|
279 | CALL TOBUF(EVTH,0)
|
---|
280 | IF (LCERFI) CALL TOBUFC(EVTH,0)
|
---|
281 | CALL COORIN(DBLE(-Z(1)))
|
---|
282 | TIME(1)=0.D0
|
---|
283 | FNPRIM =.TRUE.
|
---|
284 | IF (FPRINT) THEN
|
---|
285 | WRITE(KMPO,* )' FIRST INTERACTION AT ',EVTH(7)*0.01,' M'
|
---|
286 | END IF
|
---|
287 | END IF
|
---|
288 | CALL RMMAR(RNNO36,1,2)
|
---|
289 | GBR1=GBR11(LGLE)*GLE+GBR10(LGLE)
|
---|
290 | IF ((RNNO36.LE.GBR1).AND.(E(NP).GT.RMT2)) THEN
|
---|
291 | CALL PAIR
|
---|
292 | GO TO 982
|
---|
293 | END IF
|
---|
294 | GBR2=GBR21(LGLE)*GLE+GBR20(LGLE)
|
---|
295 | IF (RNNO36.LT.GBR2) THEN
|
---|
296 | CALL COMPT
|
---|
297 | IF((IQ(NP).NE.1))GO TO 982
|
---|
298 | GO TO1060
|
---|
299 | END IF
|
---|
300 | GBR4=GBR41(LGLE)*GLE+GBR40(LGLE)
|
---|
301 | IF (RNNO36.GE.GBR4 .AND. E(NP).GT.RMMUT2) THEN
|
---|
302 | CALL MUPAIR
|
---|
303 | GO TO 982
|
---|
304 | END IF
|
---|
305 | GBR3=GBR31(LGLE)*GLE+GBR30(LGLE)
|
---|
306 | IF (RNNO36.GE.GBR3 .AND. E(NP).GT.PITHR) THEN
|
---|
307 | CALL PIGEN
|
---|
308 | IF (NP.EQ.0) THEN
|
---|
309 | IRCODE=2
|
---|
310 | RETURN
|
---|
311 | END IF
|
---|
312 | GO TO 982
|
---|
313 | ELSE
|
---|
314 | CALL PHOTO
|
---|
315 | IF (NP.EQ.0) THEN
|
---|
316 | IRCODE=2
|
---|
317 | RETURN
|
---|
318 | END IF
|
---|
319 | IF((IQ(NP).EQ.3))GO TO 982
|
---|
320 | END IF
|
---|
321 | 1060 PEIG=E(NP)
|
---|
322 | EIG=PEIG
|
---|
323 | IF((EIG.LT.PCUT(IRL)))GO TO 970
|
---|
324 | GO TO 981
|
---|
325 | 982 CONTINUE
|
---|
326 | RETURN
|
---|
327 | 970 IF (EIG.GT.AP) THEN
|
---|
328 | IDR=1
|
---|
329 | ELSE
|
---|
330 | IDR=2
|
---|
331 | END IF
|
---|
332 | EDEP=PEIG
|
---|
333 | IRCODE=2
|
---|
334 | NP=NP-1
|
---|
335 | RETURN
|
---|
336 | 1000 EDEP=PEIG
|
---|
337 | IRCODE=2
|
---|
338 | NP=NP-1
|
---|
339 | RETURN
|
---|
340 | END
|
---|