| 1 | SUBROUTINE NKG( ENERN )
|
|---|
| 2 |
|
|---|
| 3 | C-----------------------------------------------------------------------
|
|---|
| 4 | C N(ISHIMURA) K(AMATA) G(REISEN)
|
|---|
| 5 | C
|
|---|
| 6 | C CALCULATES ELECTROMAGNETIC COMPONENT OF SHOWERS USING THE ANALYTIC
|
|---|
| 7 | C NKG FORMULAS, INCLUDING ELECTRON ENERGY THRESHOLD ELCUT(3)
|
|---|
| 8 | C SEE J.N. CAPDEVIELLE, 22ND ICRC, DUBLIN 1991, CONTRIB. HE 3.5.10
|
|---|
| 9 | C THIS SUBROUTINE IS CALLED FROM EM
|
|---|
| 10 | C ARGUMENTS:
|
|---|
| 11 | C ENERN = ENERGY OF ELECTRON/PHOTON GENERATING A SUBSHOWER
|
|---|
| 12 | C NEGATIVE FOR SUBSHOWERS TO BE SUBTRACTED AFTER
|
|---|
| 13 | C PHOTONUCLEAR REACTION
|
|---|
| 14 | C-----------------------------------------------------------------------
|
|---|
| 15 |
|
|---|
| 16 | IMPLICIT NONE
|
|---|
| 17 | *KEEP,CONST.
|
|---|
| 18 | COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
|
|---|
| 19 | DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
|
|---|
| 20 | *KEEP,ELABCT.
|
|---|
| 21 | COMMON /ELABCT/ ELCUT
|
|---|
| 22 | DOUBLE PRECISION ELCUT(4)
|
|---|
| 23 | *KEEP,NKGI.
|
|---|
| 24 | COMMON /NKGI/ SEL,SELLG,STH,ZEL,ZELLG,ZSL,DIST,
|
|---|
| 25 | * DISX,DISY,DISXY,DISYX,DLAX,DLAY,DLAXY,DLAYX,
|
|---|
| 26 | * OBSATI,RADNKG,RMOL,TLEV,TLEVCM,IALT
|
|---|
| 27 | DOUBLE PRECISION SEL(10),SELLG(10),STH(10),ZEL(10),ZELLG(10),
|
|---|
| 28 | * ZSL(10),DIST(10),
|
|---|
| 29 | * DISX(-10:10),DISY(-10:10),
|
|---|
| 30 | * DISXY(-10:10,2),DISYX(-10:10,2),
|
|---|
| 31 | * DLAX (-10:10,2),DLAY (-10:10,2),
|
|---|
| 32 | * DLAXY(-10:10,2),DLAYX(-10:10,2),
|
|---|
| 33 | * OBSATI(2),RADNKG,RMOL(2),TLEV(10),TLEVCM(10)
|
|---|
| 34 | INTEGER IALT(2)
|
|---|
| 35 | *KEEP,NKGS.
|
|---|
| 36 | COMMON /NKGS/ CZX,CZY,CZXY,CZYX,SAH,SL,ZNE
|
|---|
| 37 | DOUBLE PRECISION CZX(-10:10,2),CZY(-10:10,2),CZXY(-10:10,2),
|
|---|
| 38 | * CZYX(-10:10,2),SAH(10),SL(10),ZNE(10)
|
|---|
| 39 | *KEEP,OBSPAR.
|
|---|
| 40 | COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP,
|
|---|
| 41 | * THETPR,PHIPR,NOBSLV
|
|---|
| 42 | DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10),
|
|---|
| 43 | * THETAP,THETPR(2),PHIP,PHIPR(2)
|
|---|
| 44 | INTEGER NOBSLV
|
|---|
| 45 | *KEEP,PARPAR.
|
|---|
| 46 | COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
|
|---|
| 47 | * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
|
|---|
| 48 | DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
|
|---|
| 49 | * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
|
|---|
| 50 | INTEGER ITYPE,LEVL
|
|---|
| 51 | *KEEP,PARPAE.
|
|---|
| 52 | DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
|
|---|
| 53 | EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
|
|---|
| 54 | * (CURPAR(4), PHI ), (CURPAR(5), H ),
|
|---|
| 55 | * (CURPAR(6), T ), (CURPAR(7), X ),
|
|---|
| 56 | * (CURPAR(8), Y ), (CURPAR(9), CHI ),
|
|---|
| 57 | * (CURPAR(10),BETA), (CURPAR(11),GCM ),
|
|---|
| 58 | * (CURPAR(12),ECM )
|
|---|
| 59 | *KEEP,RUNPAR.
|
|---|
| 60 | COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
|
|---|
| 61 | * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
|
|---|
| 62 | * MONIOU,MDEBUG,NUCNUC,
|
|---|
| 63 | * CETAPE,
|
|---|
| 64 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 65 | * N1STTR,MDBASE,
|
|---|
| 66 | * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 67 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 68 | * ,GHEISH,GHESIG
|
|---|
| 69 | COMMON /RUNPAC/ DSN,HOST,USER
|
|---|
| 70 | DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
|
|---|
| 71 | REAL STEPFC
|
|---|
| 72 | INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
|
|---|
| 73 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 74 | * N1STTR,MDBASE
|
|---|
| 75 | INTEGER CETAPE
|
|---|
| 76 | CHARACTER*79 DSN
|
|---|
| 77 | CHARACTER*20 HOST,USER
|
|---|
| 78 |
|
|---|
| 79 | LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 80 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 81 | * ,GHEISH,GHESIG
|
|---|
| 82 | *KEND.
|
|---|
| 83 |
|
|---|
| 84 | DOUBLE PRECISION AE,AS,ASE,AUXIL,BS,CCP,CPC,CPCP,CPH,CSGA,
|
|---|
| 85 | * DE,DISTL,ECRI,ECR1,ECR2,ENERN,GAM,GRCUT,
|
|---|
| 86 | * G1,G2,G3,S,SC1,SC2,SIGNE,SM,SMRM,
|
|---|
| 87 | * SQRZ1I,SQZC1I,SQZC2I,SS2,SS45,TEX,THICK,THICKP,
|
|---|
| 88 | * XMOL,XNE,XS,X0,YM,YS,ZC1,ZC2,ZG1,ZG2,ZG3,Z1
|
|---|
| 89 | INTEGER IL,IOL,M
|
|---|
| 90 | EXTERNAL GAM,THICK
|
|---|
| 91 | C X0 IS RADIATON LENGTH IN AIR (G/CM**2)
|
|---|
| 92 | C (SEE ALSO MIKOCKI ET AL. J.PHYS.G.:NUCL.PART.PHYS. 17 (1991) 1303 )
|
|---|
| 93 | C GRCUT IS GREISEN CUT OFF, ECRI IS CRITICAL ENERGY IN AIR
|
|---|
| 94 | C ECR2 IS 0.4 * ECRI
|
|---|
| 95 | DATA X0 / 37.1D0 /, GRCUT / 0.1D0 /, ECRI / 0.082D0 /
|
|---|
| 96 | DATA ECR2 / 0.0328D0 /
|
|---|
| 97 | C-----------------------------------------------------------------------
|
|---|
| 98 |
|
|---|
| 99 | IF (DEBUG) WRITE(MDEBUG,*)'NKG : ',SNGL(SECPAR(1)),SNGL(ENERN)
|
|---|
| 100 |
|
|---|
| 101 | C CHECK WETHER SUBSHOWER IS SUBTRACTED
|
|---|
| 102 | IF ( ENERN .GE. 0.D0 ) THEN
|
|---|
| 103 | SIGNE = +1.D0
|
|---|
| 104 | ELSE
|
|---|
| 105 | ENERN = -ENERN
|
|---|
| 106 | SIGNE = -1.D0
|
|---|
| 107 | ENDIF
|
|---|
| 108 |
|
|---|
| 109 | C ENERGY CUT OFF IN GREISEN FORMULA
|
|---|
| 110 | C (EM PARTICLE BELOW THIS CUT CAN NOT PRODUCE A SHOWER)
|
|---|
| 111 | IF ( ENERN .LT. GRCUT ) RETURN
|
|---|
| 112 | C DON'T CALCULATE NKG FOR BACKWARD GOING PARTICLES
|
|---|
| 113 | IF ( SECPAR(3) .LE. 0.D0 ) RETURN
|
|---|
| 114 | C DON'T CALCULATE NKG IF PARTICLE BELOW THE LOWEST OBSERVATION LEVEL
|
|---|
| 115 | IF ( SECPAR(5) .LT. OBSATI(1) ) RETURN
|
|---|
| 116 |
|
|---|
| 117 | Z1 = LOG(ENERN / ECRI)
|
|---|
| 118 | SQRZ1I = 1.D0 / SQRT(Z1)
|
|---|
| 119 |
|
|---|
| 120 | C THIS CUT IS ONLY IMPORTANT FOR ELCUT > .0672
|
|---|
| 121 | ECR1 = ECR2 + ELCUT(3)
|
|---|
| 122 | IF ( ENERN .LT. ECR1 ) RETURN
|
|---|
| 123 | ZC1 = LOG(ENERN / ECR1)
|
|---|
| 124 | SQZC1I = 1.D0 / SQRT(ZC1)
|
|---|
| 125 | C LOG(ENERN/ECR2) IS LOG(ENERN / ECRI) - LOG(0.4)
|
|---|
| 126 | ZC2 = Z1 + 0.916290732D0
|
|---|
| 127 | SQZC2I = 1.D0 / SQRT(ZC2)
|
|---|
| 128 | THICKP = THICK(SECPAR(5))
|
|---|
| 129 |
|
|---|
| 130 | C LOOP OVER LEVELS
|
|---|
| 131 | DO 14 IL = 1,IALT(1)
|
|---|
| 132 | C DISREGARD LEVELS ABOVE THE PARTICLE
|
|---|
| 133 | IF ( TLEVCM(IL) .GT. SECPAR(5) ) GOTO 14
|
|---|
| 134 | C DISTANCE IN G/CM**2 .... (ALONG PHOTON-AXIS) IN RADIATION LENGTHS
|
|---|
| 135 | XMOL = (TLEV(IL) - THICKP) / ( X0 * SECPAR(3) )
|
|---|
| 136 | C CORRECT DEPTH FOR SUBSHOWERS TO BE SUBTRACTED BY 9/7
|
|---|
| 137 | IF ( SIGNE .LT. 0.D0 ) XMOL = XMOL + 1.285714286D0
|
|---|
| 138 | C XMOL IS DEPTH IN RADIATION LENGTHS
|
|---|
| 139 | IF ( XMOL .GT. 60.D0 .OR. XMOL .LT. 1.D0 ) GOTO 14
|
|---|
| 140 | C S IS AGE PARAMETER
|
|---|
| 141 | S = 3.D0 * XMOL / (XMOL + 2.D0 * Z1)
|
|---|
| 142 | IF ( S .LE. 0.2D0 ) GOTO 14
|
|---|
| 143 | SC1 = 3.D0 * XMOL / (XMOL + 2.D0 * ZC1)
|
|---|
| 144 | SC2 = 3.D0 * XMOL / (XMOL + 2.D0 * ZC2)
|
|---|
| 145 | C ELECTRON NUMBER AT OBSERVATION LEVEL
|
|---|
| 146 | CPH = .31D0 * EXP( XMOL * (1.D0 - 1.5D0 * LOG(S) ) ) * SQRZ1I
|
|---|
| 147 | CPC = EXP( XMOL * ( 1.D0 - 1.5D0 * LOG(SC1) ) ) * SQZC1I
|
|---|
| 148 | CCP = EXP( XMOL * ( 1.D0 - 1.5D0 * LOG(SC2) ) ) * SQZC2I
|
|---|
| 149 | CPCP = SIGNE * CPH * CPC / CCP
|
|---|
| 150 | C INTERMEDIATE FACTORS FOR LATERAL DISTRIBUTION AND AGE PARAMETER
|
|---|
| 151 | AE = 4.D0 * EXP( 0.915D0 * (S - 1.D0) ) / S
|
|---|
| 152 | DE = ( 1.D0 + S ) / ( 1.15D0 + 0.15D0 * S )
|
|---|
| 153 | ASE = AE**DE
|
|---|
| 154 | ZG3 = GAM( (S + 2.D0) * DE )
|
|---|
| 155 | IF ( ZG3 .LE. 0.D0 ) GOTO 14
|
|---|
| 156 | ZG1 = GAM(S * DE)
|
|---|
| 157 | ZG2 = GAM( (S + 1.D0) * DE )
|
|---|
| 158 | AUXIL = 4.D0 / (S * ASE)
|
|---|
| 159 | XNE = CPCP * ( ZG2 + AUXIL * ZG3 ) / ( ASE * (ZG1 + AUXIL*ZG2) )
|
|---|
| 160 | C SUM OF N_E AT FIXED LEVEL
|
|---|
| 161 | ZNE(IL) = ZNE(IL) + XNE
|
|---|
| 162 | SL(IL) = SL(IL) + CPCP
|
|---|
| 163 |
|
|---|
| 164 | C CALCULATE THE ELECTRON LATERAL DISTRIBUTION FOR THE 2 SELECTED
|
|---|
| 165 | C OBSERVATION LEVELS
|
|---|
| 166 | IF ( IL .EQ. IALT(1) ) THEN
|
|---|
| 167 | IOL = 1
|
|---|
| 168 | ELSEIF ( IL .EQ. IALT(2) ) THEN
|
|---|
| 169 | IOL = 2
|
|---|
| 170 | ELSE
|
|---|
| 171 | GOTO 14
|
|---|
| 172 | ENDIF
|
|---|
| 173 |
|
|---|
| 174 | C CALCULATION OF LATERAL ELECTRON DISTRIBUTION
|
|---|
| 175 | IF ( SC1 .GE. 2.25D0 ) GOTO 14
|
|---|
| 176 | G1 = GAM(4.5D0 - SC1)
|
|---|
| 177 | G2 = GAM(SC1)
|
|---|
| 178 | G3 = GAM(4.5D0 - 2.D0 * SC1)
|
|---|
| 179 | C DISTANCE IN CM BETWEEN PHOTON INITIATION AND OBSERVATION (VERTICAL)
|
|---|
| 180 | DISTL = SECPAR(5) - TLEVCM(IL)
|
|---|
| 181 | C MODULATION BY AGE PARAMETER FOLLOWING LAGUTIN & UCHAIKIN
|
|---|
| 182 | C (AGE PARAMETER LIES BETWEEN 0.2 AND 2.25)
|
|---|
| 183 | SM = 0.78D0 - 0.21D0 * SC1
|
|---|
| 184 | SMRM = 1.D0 / ( SM * RMOL(IOL) )
|
|---|
| 185 |
|
|---|
| 186 | CSGA = CPCP * SMRM**2 * G1 / ( PI2 * G2 * G3 )
|
|---|
| 187 | SS2 = SC1 - 2.D0
|
|---|
| 188 | SS45 = SC1 - 4.5D0
|
|---|
| 189 | AS = SIN( SECPAR(4) )
|
|---|
| 190 | BS = COS( SECPAR(4) )
|
|---|
| 191 | TEX = DISTL * SQRT( 1.D0 - SECPAR(3)**2 ) / SECPAR(3)
|
|---|
| 192 | C DISTANCE TO THE CENTER OF THE CASCADE (IN CM)
|
|---|
| 193 | XS = SECPAR(7) + TEX * BS - XOFF(NOBSLV+1-IOL)
|
|---|
| 194 | YS = SECPAR(8) + TEX * AS - YOFF(NOBSLV+1-IOL)
|
|---|
| 195 |
|
|---|
| 196 | C NKG-FORMULA
|
|---|
| 197 | C LOOP OVER ALL LATERAL DISTANCES GETTING THE DENSITY IN MOLIERE UNITS
|
|---|
| 198 | DO 171 M = -10,10
|
|---|
| 199 | IF ( M .EQ. 0 ) GOTO 171
|
|---|
| 200 | C X DIRECTION
|
|---|
| 201 | YM = SMRM * MAX( SQRT((DISX(M)-XS)**2 + YS**2), 1.D0 )
|
|---|
| 202 | CZX (M,IOL) = CZX (M,IOL) + CSGA * YM**SS2 * (YM+1.D0)**SS45
|
|---|
| 203 | C Y DIRECTION
|
|---|
| 204 | YM = SMRM * MAX( SQRT(XS**2 + (DISY(M)-YS)**2), 1.D0 )
|
|---|
| 205 | CZY (M,IOL) = CZY (M,IOL) + CSGA * YM**SS2 * (YM+1.D0)**SS45
|
|---|
| 206 | C XY DIRECTION
|
|---|
| 207 | YM = SMRM *
|
|---|
| 208 | * MAX( SQRT((DISXY(M,1)-XS)**2 + (DISXY(M,2)-YS)**2), 1.D0 )
|
|---|
| 209 | CZXY(M,IOL) = CZXY(M,IOL) + CSGA * YM**SS2 * (YM+1.D0)**SS45
|
|---|
| 210 | C YX DIRECTION
|
|---|
| 211 | YM = SMRM *
|
|---|
| 212 | * MAX( SQRT((DISYX(M,1)-XS)**2 + (DISYX(M,2)-YS)**2), 1.D0 )
|
|---|
| 213 | CZYX(M,IOL) = CZYX(M,IOL) + CSGA * YM**SS2 * (YM+1.D0)**SS45
|
|---|
| 214 | 171 CONTINUE
|
|---|
| 215 |
|
|---|
| 216 | 14 CONTINUE
|
|---|
| 217 |
|
|---|
| 218 | RETURN
|
|---|
| 219 | END
|
|---|
