SUBROUTINE PHOTON(IRCODE) C VERSION 4.00 -- 26 JAN 1986/1900 C****************************************************************** DOUBLE PRECISION PEIG COMMON/BOUNDS/ECUT(6),PCUT(6),VACDST *KEEP,BUFFS. COMMON /BUFFS/ RUNH,RUNE,EVTH,EVTE,DATAB,LH INTEGER MAXBUF,MAXLEN PARAMETER (MAXBUF=39*7) PARAMETER (MAXLEN=12) REAL RUNH(MAXBUF),EVTH(MAXBUF),EVTE(MAXBUF), * RUNE(MAXBUF),DATAB(MAXBUF) INTEGER LH CHARACTER*4 CRUNH,CRUNE,CEVTH,CEVTE EQUIVALENCE (RUNH(1),CRUNH), (RUNE(1),CRUNE) EQUIVALENCE (EVTH(1),CEVTH), (EVTE(1),CEVTE) *KEEP,CEREN1. COMMON /CEREN1/ CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD, * CERSIZ,LCERFI DOUBLE PRECISION CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD REAL CERSIZ LOGICAL LCERFI *KEEP,EPCONT. COMMON/EPCONT/ EDEP,RATIO,TSTEP,TUSTEP,USTEP,TVSTEP,VSTEP,IDISC, * IROLD,IRNEW,RHOFAC, EOLD,ENEW,EKE,ELKE,BETA2,GLE, * TSCAT,IAUSFL DOUBLE PRECISION EDEP,RATIO REAL TSTEP,TUSTEP,USTEP,TVSTEP,VSTEP,RHOFAC,EOLD,ENEW, * EKE,ELKE,BETA2,GLE,TSCAT INTEGER IDISC,IROLD,IRNEW,IAUSFL(29) *KEND. COMMON/GEOM/ZALTIT,BOUND(6),NEWOBS,OBSLVL(10) *KEEP,LONGI. COMMON /LONGI/ APLONG,HLONG,PLONG,SPLONG,THSTEP,THSTPI, * NSTEP,LLONGI,FLGFIT DOUBLE PRECISION APLONG(0:1040,9),HLONG(0:1024),PLONG(0:1040,9), * SPLONG(0:1040,9),THSTEP,THSTPI INTEGER NSTEP LOGICAL LLONGI,FLGFIT *KEND. COMMON /MEDIA/ NMED, RLC,RLDU,RLDUI,RHO,MSGE,MGE,MSEKE,MEKE,MLEKE, *MCMFP,MRANGE,IRAYLM,HBARO(6),HBAROI(6) CHARACTER MEDIA*24 COMMON/MEDIAC/MEDIA COMMON/MISC/KMPI,KMPO,DUNIT,NOSCAT,MED(6),RHOR(6),IRAYLR(6) DOUBLE PRECISION PRRMMU COMMON/MUON/PRRMMU,RMMU,RMMUT2 *KEEP,OBSPAR. COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP, * THETPR,PHIPR,NOBSLV DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10), * THETAP,THETPR(2),PHIP,PHIPR(2) INTEGER NOBSLV *KEEP,PARPAR. COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C, * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14), * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH INTEGER ITYPE,LEVL *KEND. COMMON/PHOTIN/EBINDA,GE0,GE1, MPGEM(1),GMFP0(500),GMFP1(500),GBR10 *(500),GBR11(500),GBR20(500),GBR21(500),GBR30(500),GBR31(500),GBR40 *(500),GBR41(500),NGR,RCO0,RCO1, RSCT0(100),RSCT1(100), COHE0(500), *COHE1(500) DOUBLE PRECISION PI0MSQ COMMON/PION/PI0MSQ,PITHR,PICMAS,PI0MAS,AMASK0,AMASKC,AMASPR,AMASNT * *KEEP,RANDPA. COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR DOUBLE PRECISION FAC,U1,U2 REAL RD(3000) INTEGER ISEED(103,10),NSEQ LOGICAL KNOR *KEEP,REJECT. COMMON /REJECT/ AVNREJ, * ALTMIN,ANEXP,THICKA,THICKD,CUTLN,EONCUT, * FNPRIM DOUBLE PRECISION AVNREJ(10) REAL ALTMIN(10),ANEXP(10),THICKA(10),THICKD(10), * CUTLN,EONCUT LOGICAL FNPRIM *KEEP,RUNPAR. COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB, * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN, * MONIOU,MDEBUG,NUCNUC, * CETAPE, * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL, * N1STTR,MDBASE, * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR, * FIX1I,FMUADD,FNKG,FPRINT,FDBASE * ,GHEISH,GHESIG COMMON /RUNPAC/ DSN,HOST,USER DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB REAL STEPFC INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC, * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL, * N1STTR,MDBASE INTEGER CETAPE CHARACTER*79 DSN CHARACTER*20 HOST,USER LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR, * FIX1I,FMUADD,FNKG,FPRINT,FDBASE * ,GHEISH,GHESIG *KEEP,STACKE. COMMON/STACKE/ E,TIME,X,Y,Z,U,V,W,DNEAR,IQ,IGEN,IR,IOBS,LPCTE,NP DOUBLE PRECISION E(60),TIME(60) REAL X(60),Y(60),Z(60),U(60),V(60),W(60),DNEAR(60) INTEGER IQ(60),IGEN(60),IR(60),IOBS(60),LPCTE(60),NP *KEND. COMMON/THRESH/RMT2,RMSQ,ESCD2,AP,API,AE,UP,UE,TE,THMOLL COMMON/UPHIOT/THETA,SINTHE,COSTHE,SINPHI, COSPHI,PI,TWOPI,PI5D2 DOUBLE PRECISION PZERO,PRM,PRMT2,RMI,VC COMMON/USEFUL/PZERO,PRM,PRMT2,RMI,VC,RM,MEDIUM,MEDOLD,IBLOBE,ICALL COMMON/ACLOCK/NCLOCK,JCLOCK DOUBLE PRECISION THICK C_____IF (NCLOCK.GT.JCLOCK) THEN C______WRITE(MDEBUG,* )' PHOTON:NP=',NP,' IR=',IR(NP),' IOBS=',IOBS(NP) C______CALL AUSGB2 C_____END IF NEWOBS=IOBS(NP) IRCODE=1 PEIG=E(NP) EIG=PEIG IRL=IR(NP) MEDIUM=MED(IRL) IF((EIG.LE.PCUT(IRL)))GO TO 970 980 CONTINUE 981 CONTINUE GLE=LOG(EIG) CALL RMMAR(RNNO35,1,2) IF ((RNNO35.EQ.0.0)) THEN RNNO35=1.E-30 END IF DPMFP=-ALOG(RNNO35) IROLD=IR(NP) 1030 CONTINUE 1031 CONTINUE IF (MEDIUM.NE.0) THEN LGLE=GE1*GLE+GE0 GMFPR0=GMFP1(LGLE)*GLE+GMFP0(LGLE) END IF 1040 CONTINUE 1041 CONTINUE IF (MEDIUM.EQ.0) THEN TSTEP=VACDST ELSE RHOFAC=RHOR(IRL)/RHO RHOFI=1./RHOFAC GMFP=GMFPR0*RHOFI IF ((IRAYLR(IRL).EQ.1)) THEN COHFAC=COHE1(LGLE)*GLE+COHE0(LGLE) GMFP=GMFP*COHFAC END IF TSTEP=GMFP*DPMFP ALTEXP=EXP(-Z(NP)*HBAROI(IRL)) TSTEP=TSTEP*ALTEXP DISC=W(NP)*TSTEP*HBAROI(IRL) IF (ABS(DISC).LT.0.065) THEN TSTEP=TSTEP*(1.-0.5*DISC*(1.-0.6666667*DISC* (1.-0.75*DISC * * (1.-0.8*DISC)))) ELSE IF(DISC.LE.-1.) THEN TSTEP=VACDST ELSE TSTEP=TSTEP/DISC*LOG(DISC+1.) END IF END IF IRNEW=IR(NP) IDISC=0 USTEP=TSTEP TUSTEP=USTEP IF((USTEP.GT.DNEAR(NP)))CALL HOWFAR IF((IDISC.GT.0))GO TO 1000 VSTEP=USTEP TVSTEP=VSTEP EDEP=PZERO USTEPU=USTEP DISC=W(NP)*USTEPU*HBAROI(IRL) USTEPU=USTEPU/ALTEXP IF (ABS(DISC).LT.0.16) THEN USTEPU=USTEPU*(1.+.5*DISC*(1.+.33333333*DISC* (1.+0.25*DISC* ( * 1.+0.2*DISC)))) ELSE USTEPU=USTEPU/DISC*(EXP(DISC)-1.) END IF X(NP)=X(NP)+U(NP)*USTEP Y(NP)=Y(NP)+V(NP)*USTEP ZOLD =Z(NP) Z(NP)=Z(NP)+W(NP)*USTEP TIME(NP)=TIME(NP)+TVSTEP*VC C ADD PHOTONS TO THE LONGITUDINAL DEVELOPMENT IF ( LLONGI ) THEN C FIND FIRST THE EQUIVALENT LEVELS C IF STARTING POINT BELOW LOWEST LEVEL THEN DON'T CHECK IF ( HLONG(NSTEP) .LE. -ZOLD ) THEN LPCT1 = LPCTE(NP) C Z NEW IS PROBABLY ONLY LITTLE BELOW Z OLD, THEREFORE INCREMENTAL SEARCH DO 6002 I1 = LPCT1,NSTEP IF ( HLONG(I1) .LT. -Z(NP) ) GOTO 6003 6002 CONTINUE I1 = NSTEP + 1 6003 CONTINUE LPCT2 = I1 - 1 DO 485 I=LPCT1,LPCT2 PLONG(I,1) = PLONG(I,1) + 1.D0 485 CONTINUE LPCTE(NP) = LPCT2 + 1 ENDIF ENDIF DNEAR(NP)=DNEAR(NP)-USTEP IF (MEDIUM.NE.0) THEN DPMFP=MAX(0.,DPMFP-USTEPU/GMFP) END IF IROLD=IR(NP) MEDOLD=MEDIUM IF (IRNEW.NE.IROLD) THEN IR(NP)=IRNEW IRL=IRNEW MEDIUM=MED(IRL) IF((EIG.LE.PCUT(IRL)))GO TO 970 END IF IF (NEWOBS.GT.IOBS(NP)) THEN CALL AUSGAB IOBS(NP)=NEWOBS END IF IF((IDISC.LT.0))GO TO 1000 IF((MEDIUM.NE.MEDOLD))GO TO 1042 IF((MEDIUM.NE.0.AND.DPMFP.LE.1.E-6))GO TO 1032 GO TO 1041 1042 CONTINUE GO TO 1031 1032 CONTINUE IF ((IRAYLR(IRL).EQ.1)) THEN CALL RMMAR(RNNO37,1,2) IF ((RNNO37.LE.(1.0-COHFAC))) THEN 1050 CONTINUE 1051 CONTINUE CALL RMMAR(XXX,1,2) LXXX=RCO1*XXX+RCO0 X2=RSCT1(LXXX)*XXX+RSCT0(LXXX) Q2=X2*RMSQ*.23547885E-02 COSTHE=1.-Q2/(2.*E(NP)*E(NP)) IF((ABS(COSTHE).GT.1.0))GO TO 1050 CSQTHE=COSTHE*COSTHE REJF=(1.0+CSQTHE)*.5 CALL RMMAR(RNNORJ,1,2) IF((RNNORJ.LE.REJF))GO TO1052 GO TO 1051 1052 CONTINUE SINTHE=SQRT(AMAX1(0.,1.0-CSQTHE)) CALL UPHI(2,1) GOTO 980 END IF END IF IF ( .NOT. FNPRIM ) THEN X(1)=0. Y(1)=0. EVTH(5)=X(1) EVTH(6)=-Y(1) IF (FIX1I) THEN Z(1)=-FIXHEI NP=1 LPCTE(1)=MIN(NSTEP,INT(THICK(FIXHEI)*THSTPI)+1) SITHET=SQRT(1.D0-SECPAR(3)**2) U(1)=SITHET*COS(-SECPAR(4)) V(1)=SITHET*SIN(-SECPAR(4)) W(1)=SECPAR(3) RADINV=1.5-0.5*(U(1)**2+V(1)**2+W(1)**2) U(1)=U(1)*RADINV V(1)=V(1)*RADINV W(1)=W(1)*RADINV END IF EVTH(7)=-Z(1) CALL TOBUF(EVTH,0) IF (LCERFI) CALL TOBUFC(EVTH,0) CALL COORIN(DBLE(-Z(1))) TIME(1)=0.D0 FNPRIM =.TRUE. IF (FPRINT) THEN WRITE(KMPO,* )' FIRST INTERACTION AT ',EVTH(7)*0.01,' M' END IF END IF CALL RMMAR(RNNO36,1,2) GBR1=GBR11(LGLE)*GLE+GBR10(LGLE) IF ((RNNO36.LE.GBR1).AND.(E(NP).GT.RMT2)) THEN CALL PAIR GO TO 982 END IF GBR2=GBR21(LGLE)*GLE+GBR20(LGLE) IF (RNNO36.LT.GBR2) THEN CALL COMPT IF((IQ(NP).NE.1))GO TO 982 GO TO1060 END IF GBR4=GBR41(LGLE)*GLE+GBR40(LGLE) IF (RNNO36.GE.GBR4 .AND. E(NP).GT.RMMUT2) THEN CALL MUPAIR GO TO 982 END IF GBR3=GBR31(LGLE)*GLE+GBR30(LGLE) IF (RNNO36.GE.GBR3 .AND. E(NP).GT.PITHR) THEN CALL PIGEN IF (NP.EQ.0) THEN IRCODE=2 RETURN END IF GO TO 982 ELSE CALL PHOTO IF (NP.EQ.0) THEN IRCODE=2 RETURN END IF IF((IQ(NP).EQ.3))GO TO 982 END IF 1060 PEIG=E(NP) EIG=PEIG IF((EIG.LT.PCUT(IRL)))GO TO 970 GO TO 981 982 CONTINUE RETURN 970 IF (EIG.GT.AP) THEN IDR=1 ELSE IDR=2 END IF EDEP=PEIG IRCODE=2 NP=NP-1 RETURN 1000 EDEP=PEIG IRCODE=2 NP=NP-1 RETURN END