| 1 | SUBROUTINE EGSINI
|
|---|
| 2 |
|
|---|
| 3 | C---------------------------------------------------------------------
|
|---|
| 4 | C E(LECTRON) G(AMMA) S(HOWER) INI(TIALIZATION)
|
|---|
| 5 | C
|
|---|
| 6 | C INITIALIZES EGS4 PACKAGE
|
|---|
| 7 | C THIS SUBROUTINE IS CALLED FROM START
|
|---|
| 8 | C
|
|---|
| 9 | C DESIGN : D. HECK IK3 FZK KARLSRUHE
|
|---|
| 10 | C---------------------------------------------------------------------
|
|---|
| 11 |
|
|---|
| 12 | *KEEP,ATMOS.
|
|---|
| 13 | COMMON /ATMOS/ AATM,BATM,CATM,DATM
|
|---|
| 14 | DOUBLE PRECISION AATM(5),BATM(5),CATM(5),DATM(5)
|
|---|
| 15 | *KEND.
|
|---|
| 16 | COMMON/BOUNDS/ECUT(6),PCUT(6),VACDST
|
|---|
| 17 | *KEEP,ELABCT.
|
|---|
| 18 | COMMON /ELABCT/ ELCUT
|
|---|
| 19 | DOUBLE PRECISION ELCUT(4)
|
|---|
| 20 | *KEND.
|
|---|
| 21 | COMMON/GEOM/ZALTIT,BOUND(6),NEWOBS,OBSLVL(10)
|
|---|
| 22 | COMMON /MEDIA/ NMED, RLC,RLDU,RLDUI,RHO,MSGE,MGE,MSEKE,MEKE,MLEKE,
|
|---|
| 23 | *MCMFP,MRANGE,IRAYLM,HBARO(6),HBAROI(6)
|
|---|
| 24 | CHARACTER MEDIA*24
|
|---|
| 25 | COMMON/MEDIAC/MEDIA
|
|---|
| 26 | DOUBLE PRECISION PRRMMU
|
|---|
| 27 | COMMON/MUON/PRRMMU,RMMU,RMMUT2
|
|---|
| 28 | COMMON/MISC/KMPI,KMPO,DUNIT,NOSCAT,MED(6),RHOR(6),IRAYLR(6)
|
|---|
| 29 | *KEEP,OBSPAR.
|
|---|
| 30 | COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP,
|
|---|
| 31 | * THETPR,PHIPR,NOBSLV
|
|---|
| 32 | DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10),
|
|---|
| 33 | * THETAP,THETPR(2),PHIP,PHIPR(2)
|
|---|
| 34 | INTEGER NOBSLV
|
|---|
| 35 | *KEEP,PAM.
|
|---|
| 36 | COMMON /PAM/ PAMA,SIGNUM
|
|---|
| 37 | DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
|
|---|
| 38 | *KEEP,PARPAR.
|
|---|
| 39 | COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
|
|---|
| 40 | * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
|
|---|
| 41 | DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
|
|---|
| 42 | * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
|
|---|
| 43 | INTEGER ITYPE,LEVL
|
|---|
| 44 | *KEND.
|
|---|
| 45 | DOUBLE PRECISION PI0MSQ
|
|---|
| 46 | COMMON/PION/PI0MSQ,PITHR,PICMAS,PI0MAS,AMASK0,AMASKC,AMASPR,AMASNT
|
|---|
| 47 | *KEEP,REJECT.
|
|---|
| 48 | COMMON /REJECT/ AVNREJ,
|
|---|
| 49 | * ALTMIN,ANEXP,THICKA,THICKD,CUTLN,EONCUT,
|
|---|
| 50 | * FNPRIM
|
|---|
| 51 | DOUBLE PRECISION AVNREJ(10)
|
|---|
| 52 | REAL ALTMIN(10),ANEXP(10),THICKA(10),THICKD(10),
|
|---|
| 53 | * CUTLN,EONCUT
|
|---|
| 54 | LOGICAL FNPRIM
|
|---|
| 55 | *KEEP,RUNPAR.
|
|---|
| 56 | COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
|
|---|
| 57 | * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
|
|---|
| 58 | * MONIOU,MDEBUG,NUCNUC,
|
|---|
| 59 | * CETAPE,
|
|---|
| 60 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 61 | * N1STTR,MDBASE,
|
|---|
| 62 | * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 63 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 64 | * ,GHEISH,GHESIG
|
|---|
| 65 | COMMON /RUNPAC/ DSN,HOST,USER
|
|---|
| 66 | DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
|
|---|
| 67 | REAL STEPFC
|
|---|
| 68 | INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
|
|---|
| 69 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 70 | * N1STTR,MDBASE
|
|---|
| 71 | INTEGER CETAPE
|
|---|
| 72 | CHARACTER*79 DSN
|
|---|
| 73 | CHARACTER*20 HOST,USER
|
|---|
| 74 |
|
|---|
| 75 | LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 76 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 77 | * ,GHEISH,GHESIG
|
|---|
| 78 | *KEND.
|
|---|
| 79 | COMMON/THRESH/RMT2,RMSQ,ESCD2,AP,API,AE,UP,UE,TE,THMOLL
|
|---|
| 80 | DOUBLE PRECISION PZERO,PRM,PRMT2,RMI,VC
|
|---|
| 81 | COMMON/USEFUL/PZERO,PRM,PRMT2,RMI,VC,RM,MEDIUM,MEDOLD,IBLOBE,ICALL
|
|---|
| 82 | COMMON/ACLOCK/NCLOCK,JCLOCK
|
|---|
| 83 | CHARACTER MEDARR*24
|
|---|
| 84 | DATA MEDARR/'AIR-NTP '/
|
|---|
| 85 | C---------------------------------------------------------------------
|
|---|
| 86 |
|
|---|
| 87 | IF((DEBUG))WRITE(MDEBUG,* )'EGSINI:'
|
|---|
| 88 | C*** INITIALISATION BEFORE THE FIRST CALL OF EGS4
|
|---|
| 89 | ICALL = 1
|
|---|
| 90 | IF (( DEBUG )) THEN
|
|---|
| 91 | KMPO = MDEBUG
|
|---|
| 92 | ELSE
|
|---|
| 93 | KMPO = MONIOU
|
|---|
| 94 | END IF
|
|---|
| 95 | WRITE(KMPO,10)
|
|---|
| 96 | 10 FORMAT (' START AIR SHOWER SUBROUTINE EGS4 VERSION 64S (OCT 94)'/
|
|---|
| 97 | * ' MU+MU- PAIRS ARE FORMED IN ANALOGY WITH E+E- PAIRS. THIS IS ',
|
|---|
| 98 | * 'CERTAINLY NOT'/' CORRECT BECAUSE OF DIFFERENT FORM FACTORS. T',
|
|---|
| 99 | * 'HEREFORE TREAT MUON PAIRS WITH CAUTION.')
|
|---|
| 100 | C*** SET PARTICLE MASSES AND PHYSICAL CONSTANTS
|
|---|
| 101 | RM = PAMA(2)*1.D3
|
|---|
| 102 | RMT2 = 2.*RM
|
|---|
| 103 | RMSQ = RM**2
|
|---|
| 104 | PRRMMU = PAMA(5)*1.D3
|
|---|
| 105 | RMMU = PRRMMU
|
|---|
| 106 | RMMUT2 = 2.*RMMU
|
|---|
| 107 | PICMAS = PAMA(8)*1.D3
|
|---|
| 108 | PI0MAS = PAMA(7)*1.D3
|
|---|
| 109 | PI0MSQ = PI0MAS**2
|
|---|
| 110 | AMASKC = PAMA(11)*1.D3
|
|---|
| 111 | AMASK0 = PAMA(10)*1.D3
|
|---|
| 112 | AMASPR = PAMA(14)*1.D3
|
|---|
| 113 | AMASNT = PAMA(13)*1.D3
|
|---|
| 114 | C*** INVERSE OF VELOCITY OF LIGHT
|
|---|
| 115 | VC = 1.D0/C(25)
|
|---|
| 116 | C*** PION-PRODUCTION THRESHOLD
|
|---|
| 117 | PITHR = 152.0
|
|---|
| 118 | C*** NMED AND DUNIT DEFAULT TO 1,I.E. ONE MEDIUM AND WE WORK IN CM
|
|---|
| 119 | MEDIUM=1
|
|---|
| 120 | DO 21 I=1,24
|
|---|
| 121 | MEDIA(I:I)=MEDARR(I:I)
|
|---|
| 122 | 21 CONTINUE
|
|---|
| 123 | 22 CONTINUE
|
|---|
| 124 | C*** BOUNDARY 1= TOP OF ATMOSPHERE (SEE SUBROUTINE HOWFAR)
|
|---|
| 125 | BOUND(1)=11300000.
|
|---|
| 126 | BOUND(2)= 4000000.
|
|---|
| 127 | BOUND(3)= 1000000.
|
|---|
| 128 | BOUND(4)= 400000.
|
|---|
| 129 | BOUND(5)= 0.
|
|---|
| 130 | BOUND(6)= -0.
|
|---|
| 131 | MED(1)=0
|
|---|
| 132 | MED(6)=0
|
|---|
| 133 | C*** VACUUM IN REGIONS 1 AND 6, AIR IN REGION 2 TO 5
|
|---|
| 134 | DO 31 IRL=2,5
|
|---|
| 135 | MED(IRL)= 1
|
|---|
| 136 | C *** VALUE OF BARO-HEIGHT FROM HILLAS, PRIV. COMMUN. (AUG 1988)
|
|---|
| 137 | HBARO(IRL)=CATM(6-IRL)
|
|---|
| 138 | HBAROI(IRL)=1./HBARO(IRL)
|
|---|
| 139 | RHOR(IRL) =BATM(6-IRL)*HBAROI(IRL)
|
|---|
| 140 | C *** NEEDED FOR REGION 2 TO 5 SINCE NO TRANSPORT ELSEWHERE
|
|---|
| 141 | C *** ECUT IS TOTAL ENERGY
|
|---|
| 142 | C *** TERMINATE ELECTRON HISTORIES AT ECUT (GEV TO MEV CONVERTED)
|
|---|
| 143 | ECUT(IRL)=1000.D0*ELCUT(3)+RM
|
|---|
| 144 | C *** TERMINATE PHOTON HISTORIES AT PCUT (GEV TO MEV CONVERTED)
|
|---|
| 145 | PCUT(IRL)=1000.D0*ELCUT(4)
|
|---|
| 146 | 31 CONTINUE
|
|---|
| 147 | 32 CONTINUE
|
|---|
| 148 | OPEN(UNIT=KMPI,FILE='EGSDAT2',STATUS='OLD')
|
|---|
| 149 | C*** PICK UP CROSS SECTION DATA FOR AIR-NTP FROM UNIT KMPI=12
|
|---|
| 150 | CALL HATCH
|
|---|
| 151 | CLOSE(UNIT=KMPI)
|
|---|
| 152 | C*** INVERTED PHOTON THRESHOLD
|
|---|
| 153 | API=1./AP
|
|---|
| 154 | WRITE(KMPO,40) (AE-RM)*.001,AP*.001,ECUT(2)*.001,PCUT(2)*.001
|
|---|
| 155 | 40 FORMAT (' ELECTRONS CAN BE CREATED AND ANY ELECTRON FOLLOWED DO',
|
|---|
| 156 | * 'WN TO'/T40,F8.4,' GEV KINETIC ENERGY'/' GAMMAS CAN BE CREATED',
|
|---|
| 157 | * ' AND ANY GAMMA FOLLOWED DOWN TO'/T40,F8.4,' GEV ENERGY'/' ELE',
|
|---|
| 158 | * 'CTRON HISTORIES ARE TERMINATED AT',F10.4,' GEV'/' GAMMA HISTO',
|
|---|
| 159 | * 'RIES ARE TERMINATED AT ',F10.4,' GEV'/)
|
|---|
| 160 | cc IF((DEBUG.OR.(JCLOCK.GT.1)))WRITE(KMPO,50)
|
|---|
| 161 | cc50 FORMAT (7X,' PART|TOT.ENERGY|ANGLE Z|ANGLE X|ALTITUDE|',
|
|---|
| 162 | cc * ' TIME | POS. X | POS. Y |GENER|',/,8X,'ICLE|',
|
|---|
| 163 | cc * ' (GEV) |COSTHET| (RAD) | (CM) | (MSEC) | (CM) |',
|
|---|
| 164 | cc * ' (CM) |ATION|')
|
|---|
| 165 | C*** CALCULATE THE LAYER THICKNESS BELOW EACH DETECTOR
|
|---|
| 166 | DO 61 IDET=1,NOBSLV
|
|---|
| 167 | C *** NECESSARY BECAUSE OF DOUBLE PRECIS.
|
|---|
| 168 | OBSLVL(IDET)=OBSLEV(IDET)
|
|---|
| 169 | DO 71 JREG=2,5
|
|---|
| 170 | IF (OBSLVL(IDET).GE.BOUND(JREG)) THEN
|
|---|
| 171 | KREG=JREG
|
|---|
| 172 | GO TO 80
|
|---|
| 173 | END IF
|
|---|
| 174 | 71 CONTINUE
|
|---|
| 175 | 72 CONTINUE
|
|---|
| 176 | WRITE(KMPO,90) IDET,OBSLVL(IDET)*0.01
|
|---|
| 177 | 90 FORMAT (' EGS4 :', ' DETECTOR ',I2,' AT ',E10.3,' M IS OUT OF A
|
|---|
| 178 | *TMOSPHERE')
|
|---|
| 179 | STOP
|
|---|
| 180 | 80 THICKD(IDET)=EXP(-OBSLVL(IDET)*HBAROI(KREG))
|
|---|
| 181 | THICKA(IDET)=RHOR(KREG)*HBARO(KREG)*(1.-THICKD(IDET))
|
|---|
| 182 | C *** MIN ALTITUDE FOR REJECT IS OBSERVATION LEVEL+3*36.6 G/CM**2
|
|---|
| 183 | ALTMIN(IDET)=-HBARO(KREG)*LOG(MAX(1.E-37,(1.-(THICKA(IDET)+109.8
|
|---|
| 184 | * )* HBAROI(KREG)/RHOR(KREG))))
|
|---|
| 185 | ALTMIN(IDET)=MIN(ALTMIN(IDET),BOUND(1))
|
|---|
| 186 | 61 CONTINUE
|
|---|
| 187 | 62 CONTINUE
|
|---|
| 188 | RETURN
|
|---|
| 189 | END
|
|---|
