| 1 |
|
|---|
| 2 | PROGRAM MAIN
|
|---|
| 3 |
|
|---|
| 4 | C-----------------------------------------------------------------------
|
|---|
| 5 | C MAIN PROGRAM
|
|---|
| 6 | C
|
|---|
| 7 | C SIMULATION OF EXTENSIVE AIR SHOWERS
|
|---|
| 8 | C PREPARES INITIALIZATIONS
|
|---|
| 9 | C GENERATES SHOWERS IN THE SHOWER LOOP
|
|---|
| 10 | C TREATES PARTICLES IN THE PARTICLE LOOP
|
|---|
| 11 | C PERFORMS PRINTING OF TABLES AT END OF SHOWER AND AT END OF RUN
|
|---|
| 12 | C-----------------------------------------------------------------------
|
|---|
| 13 |
|
|---|
| 14 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 15 | parameter (xct=1)
|
|---|
| 16 | parameter (yct=2)
|
|---|
| 17 | parameter (zct=3)
|
|---|
| 18 | parameter (ctthet=4)
|
|---|
| 19 | parameter (ctphi=5)
|
|---|
| 20 | parameter (ctdiam=6)
|
|---|
| 21 | parameter (ctfoc=7)
|
|---|
| 22 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 23 |
|
|---|
| 24 | IMPLICIT DOUBLE PRECISION (A-H,O-Z)
|
|---|
| 25 | *KEEP,BAL.
|
|---|
| 26 | COMMON /BAL/ EBAL
|
|---|
| 27 | DOUBLE PRECISION EBAL(10)
|
|---|
| 28 | *KEEP,BUFFS.
|
|---|
| 29 | COMMON /BUFFS/ RUNH,RUNE,EVTH,EVTE,DATAB,LH
|
|---|
| 30 | INTEGER MAXBUF,MAXLEN
|
|---|
| 31 | PARAMETER (MAXBUF=39*7)
|
|---|
| 32 | PARAMETER (MAXLEN=12)
|
|---|
| 33 | REAL RUNH(MAXBUF),EVTH(MAXBUF),EVTE(MAXBUF),
|
|---|
| 34 | * RUNE(MAXBUF),DATAB(MAXBUF)
|
|---|
| 35 | INTEGER LH
|
|---|
| 36 | CHARACTER*4 CRUNH,CRUNE,CEVTH,CEVTE
|
|---|
| 37 | EQUIVALENCE (RUNH(1),CRUNH), (RUNE(1),CRUNE)
|
|---|
| 38 | EQUIVALENCE (EVTH(1),CEVTH), (EVTE(1),CEVTE)
|
|---|
| 39 | *KEEP,CHISTA.
|
|---|
| 40 | COMMON /CHISTA/ IHYCHI,IKACHI,IMUCHI,INNCHI,INUCHI,IPICHI
|
|---|
| 41 | INTEGER IHYCHI(124),IKACHI(124),IMUCHI(124),
|
|---|
| 42 | * INNCHI(124),INUCHI(124),IPICHI(124)
|
|---|
| 43 | *KEEP,CONST.
|
|---|
| 44 | COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
|
|---|
| 45 | DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
|
|---|
| 46 | *KEEP,CURVE.
|
|---|
| 47 | COMMON /CURVE/ CHAPAR,DEP,ERR,NSTP
|
|---|
| 48 | DOUBLE PRECISION CHAPAR(1100),DEP(1100),ERR(1100)
|
|---|
| 49 | INTEGER NSTP
|
|---|
| 50 | *KEEP,ELADPM.
|
|---|
| 51 | COMMON /ELADPM/ ELMEAN,ELMEAA,IELDPM,IELDPA
|
|---|
| 52 | DOUBLE PRECISION ELMEAN(37),ELMEAA(37)
|
|---|
| 53 | INTEGER IELDPM(37,13),IELDPA(37,13)
|
|---|
| 54 | *KEEP,ELASTY.
|
|---|
| 55 | COMMON /ELASTY/ ELAST,IELIS,IELHM,IELNU,IELPI
|
|---|
| 56 | DOUBLE PRECISION ELAST
|
|---|
| 57 | INTEGER IELIS(20),IELHM(20),IELNU(20),IELPI(20)
|
|---|
| 58 | *KEEP,GENER.
|
|---|
| 59 | COMMON /GENER/ GEN,ALEVEL
|
|---|
| 60 | DOUBLE PRECISION GEN,ALEVEL
|
|---|
| 61 | *KEEP,IRET.
|
|---|
| 62 | COMMON /IRET/ IRET1,IRET2
|
|---|
| 63 | INTEGER IRET1,IRET2
|
|---|
| 64 | *KEEP,ISTA.
|
|---|
| 65 | COMMON /ISTA/ IFINET,IFINNU,IFINKA,IFINPI,IFINHY
|
|---|
| 66 | INTEGER IFINET,IFINNU,IFINKA,IFINPI,IFINHY
|
|---|
| 67 | *KEEP,LONGI.
|
|---|
| 68 | COMMON /LONGI/ APLONG,HLONG,PLONG,SPLONG,THSTEP,THSTPI,
|
|---|
| 69 | * NSTEP,LLONGI,FLGFIT
|
|---|
| 70 | DOUBLE PRECISION APLONG(0:1040,9),HLONG(0:1024),PLONG(0:1040,9),
|
|---|
| 71 | * SPLONG(0:1040,9),THSTEP,THSTPI
|
|---|
| 72 | INTEGER NSTEP
|
|---|
| 73 | LOGICAL LLONGI,FLGFIT
|
|---|
| 74 | *KEEP,MPARTI.
|
|---|
| 75 | COMMON /MPARTI/ MPARTO
|
|---|
| 76 | DOUBLE PRECISION MPARTO(10,25),MPHOTO(10),MPOSIT(10),MELECT(10),
|
|---|
| 77 | * MNU(10),MMUP(10),MMUM(10),MPI0(10),MPIP(10),
|
|---|
| 78 | * MPIM(10),MK0L(10),MKPL(10),MKMI(10),MNEUTR(10),
|
|---|
| 79 | * MPROTO(10),MPROTB(10),MK0S(10),MHYP(10),
|
|---|
| 80 | * MNEUTB(10),MDEUT(10),MTRIT(10),MALPHA(10),
|
|---|
| 81 | * MOTHER(10)
|
|---|
| 82 | EQUIVALENCE (MPARTO(1, 1),MPHOTO(1)), (MPARTO(1, 2),MPOSIT(1)),
|
|---|
| 83 | * (MPARTO(1, 3),MELECT(1)), (MPARTO(1, 4),MNU(1)) ,
|
|---|
| 84 | * (MPARTO(1, 5),MMUP(1)) , (MPARTO(1, 6),MMUM(1)) ,
|
|---|
| 85 | * (MPARTO(1, 7),MPI0(1)) , (MPARTO(1, 8),MPIP(1)) ,
|
|---|
| 86 | * (MPARTO(1, 9),MPIM(1)) , (MPARTO(1,10),MK0L(1)) ,
|
|---|
| 87 | * (MPARTO(1,11),MKPL(1)) , (MPARTO(1,12),MKMI(1)) ,
|
|---|
| 88 | * (MPARTO(1,13),MNEUTR(1)), (MPARTO(1,14),MPROTO(1)),
|
|---|
| 89 | * (MPARTO(1,15),MPROTB(1)), (MPARTO(1,16),MK0S(1)) ,
|
|---|
| 90 | * (MPARTO(1,18),MHYP(1)) , (MPARTO(1,19),MDEUT(1)) ,
|
|---|
| 91 | * (MPARTO(1,20),MTRIT(1)) , (MPARTO(1,21),MALPHA(1)),
|
|---|
| 92 | * (MPARTO(1,22),MOTHER(1)), (MPARTO(1,25),MNEUTB(1))
|
|---|
| 93 | *KEEP,MULT.
|
|---|
| 94 | COMMON /MULT/ EKINL,MSMM,MULTMA,MULTOT
|
|---|
| 95 | DOUBLE PRECISION EKINL
|
|---|
| 96 | INTEGER MSMM,MULTMA(37,13),MULTOT(37,13)
|
|---|
| 97 | *KEEP,MUPART.
|
|---|
| 98 | COMMON /MUPART/ AMUPAR,BCUT,CMUON,FMUBRM,FMUORG
|
|---|
| 99 | DOUBLE PRECISION AMUPAR(14),BCUT,CMUON(11)
|
|---|
| 100 | LOGICAL FMUBRM,FMUORG
|
|---|
| 101 | *KEEP,NCOUNT.
|
|---|
| 102 | COMMON /NCOUNT/ NCOUN
|
|---|
| 103 | INTEGER NCOUN(8)
|
|---|
| 104 | *KEEP,NKGI.
|
|---|
| 105 | COMMON /NKGI/ SEL,SELLG,STH,ZEL,ZELLG,ZSL,DIST,
|
|---|
| 106 | * DISX,DISY,DISXY,DISYX,DLAX,DLAY,DLAXY,DLAYX,
|
|---|
| 107 | * OBSATI,RADNKG,RMOL,TLEV,TLEVCM,IALT
|
|---|
| 108 | DOUBLE PRECISION SEL(10),SELLG(10),STH(10),ZEL(10),ZELLG(10),
|
|---|
| 109 | * ZSL(10),DIST(10),
|
|---|
| 110 | * DISX(-10:10),DISY(-10:10),
|
|---|
| 111 | * DISXY(-10:10,2),DISYX(-10:10,2),
|
|---|
| 112 | * DLAX (-10:10,2),DLAY (-10:10,2),
|
|---|
| 113 | * DLAXY(-10:10,2),DLAYX(-10:10,2),
|
|---|
| 114 | * OBSATI(2),RADNKG,RMOL(2),TLEV(10),TLEVCM(10)
|
|---|
| 115 | INTEGER IALT(2)
|
|---|
| 116 | *KEEP,NKGS.
|
|---|
| 117 | COMMON /NKGS/ CZX,CZY,CZXY,CZYX,SAH,SL,ZNE
|
|---|
| 118 | DOUBLE PRECISION CZX(-10:10,2),CZY(-10:10,2),CZXY(-10:10,2),
|
|---|
| 119 | * CZYX(-10:10,2),SAH(10),SL(10),ZNE(10)
|
|---|
| 120 | *KEEP,NPARTI.
|
|---|
| 121 | COMMON /NPARTI/ NPARTO,MUOND
|
|---|
| 122 | DOUBLE PRECISION NPARTO(10,25),NPHOTO(10),NPOSIT(10),NELECT(10),
|
|---|
| 123 | * NNU(10),NMUP(10),NMUM(10),NPI0(10),NPIP(10),
|
|---|
| 124 | * NPIM(10),NK0L(10),NKPL(10),NKMI(10),NNEUTR(10),
|
|---|
| 125 | * NPROTO(10),NPROTB(10),NK0S(10),NHYP(10),
|
|---|
| 126 | * NNEUTB(10),NDEUT(10),NTRIT(10),NALPHA(10),
|
|---|
| 127 | * NOTHER(10),MUOND
|
|---|
| 128 | EQUIVALENCE (NPARTO(1, 1),NPHOTO(1)), (NPARTO(1, 2),NPOSIT(1)),
|
|---|
| 129 | * (NPARTO(1, 3),NELECT(1)), (NPARTO(1, 4),NNU(1)) ,
|
|---|
| 130 | * (NPARTO(1, 5),NMUP(1)) , (NPARTO(1, 6),NMUM(1)) ,
|
|---|
| 131 | * (NPARTO(1, 7),NPI0(1)) , (NPARTO(1, 8),NPIP(1)) ,
|
|---|
| 132 | * (NPARTO(1, 9),NPIM(1)) , (NPARTO(1,10),NK0L(1)) ,
|
|---|
| 133 | * (NPARTO(1,11),NKPL(1)) , (NPARTO(1,12),NKMI(1)) ,
|
|---|
| 134 | * (NPARTO(1,13),NNEUTR(1)), (NPARTO(1,14),NPROTO(1)),
|
|---|
| 135 | * (NPARTO(1,15),NPROTB(1)), (NPARTO(1,16),NK0S(1)) ,
|
|---|
| 136 | * (NPARTO(1,18),NHYP(1)) , (NPARTO(1,19),NDEUT(1)) ,
|
|---|
| 137 | * (NPARTO(1,20),NTRIT(1)) , (NPARTO(1,21),NALPHA(1)),
|
|---|
| 138 | * (NPARTO(1,22),NOTHER(1)), (NPARTO(1,25),NNEUTB(1))
|
|---|
| 139 | *KEEP,OBSPAR.
|
|---|
| 140 | COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP,
|
|---|
| 141 | * THETPR,PHIPR,NOBSLV
|
|---|
| 142 | DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10),
|
|---|
| 143 | * THETAP,THETPR(2),PHIP,PHIPR(2)
|
|---|
| 144 | INTEGER NOBSLV
|
|---|
| 145 | *KEEP,PAM.
|
|---|
| 146 | COMMON /PAM/ PAMA,SIGNUM
|
|---|
| 147 | DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
|
|---|
| 148 | *KEEP,PARPAR.
|
|---|
| 149 | COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
|
|---|
| 150 | * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
|
|---|
| 151 | DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
|
|---|
| 152 | * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
|
|---|
| 153 | INTEGER ITYPE,LEVL
|
|---|
| 154 | *KEEP,PARPAE.
|
|---|
| 155 | DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
|
|---|
| 156 | EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
|
|---|
| 157 | * (CURPAR(4), PHI ), (CURPAR(5), H ),
|
|---|
| 158 | * (CURPAR(6), T ), (CURPAR(7), X ),
|
|---|
| 159 | * (CURPAR(8), Y ), (CURPAR(9), CHI ),
|
|---|
| 160 | * (CURPAR(10),BETA), (CURPAR(11),GCM ),
|
|---|
| 161 | * (CURPAR(12),ECM )
|
|---|
| 162 | *KEEP,PBALA.
|
|---|
| 163 | COMMON /PBALA/ PBAL
|
|---|
| 164 | DOUBLE PRECISION PBAL(10)
|
|---|
| 165 | *KEEP,PRIMSP.
|
|---|
| 166 | COMMON /PRIMSP/ PSLOPE,LLIMIT,ULIMIT,LL,UL,SLEX,ISPEC
|
|---|
| 167 | DOUBLE PRECISION PSLOPE,LLIMIT,ULIMIT,LL,UL,SLEX
|
|---|
| 168 | INTEGER ISPEC
|
|---|
| 169 | *KEEP,RANDPA.
|
|---|
| 170 | COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
|
|---|
| 171 | DOUBLE PRECISION FAC,U1,U2
|
|---|
| 172 | REAL RD(3000)
|
|---|
| 173 | INTEGER ISEED(103,10),NSEQ
|
|---|
| 174 | LOGICAL KNOR
|
|---|
| 175 | *KEEP,RECORD.
|
|---|
| 176 | COMMON /RECORD/ IRECOR
|
|---|
| 177 | INTEGER IRECOR
|
|---|
| 178 | *KEEP,REJECT.
|
|---|
| 179 | COMMON /REJECT/ AVNREJ,
|
|---|
| 180 | * ALTMIN,ANEXP,THICKA,THICKD,CUTLN,EONCUT,
|
|---|
| 181 | * FNPRIM
|
|---|
| 182 | DOUBLE PRECISION AVNREJ(10)
|
|---|
| 183 | REAL ALTMIN(10),ANEXP(10),THICKA(10),THICKD(10),
|
|---|
| 184 | * CUTLN,EONCUT
|
|---|
| 185 | LOGICAL FNPRIM
|
|---|
| 186 | *KEEP,RESON.
|
|---|
| 187 | COMMON /RESON/ RDRES,RESRAN,IRESPAR
|
|---|
| 188 | REAL RDRES(2),RESRAN(1000)
|
|---|
| 189 | INTEGER IRESPAR
|
|---|
| 190 |
|
|---|
| 191 | *KEEP,RUNPAR.
|
|---|
| 192 | COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
|
|---|
| 193 | * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
|
|---|
| 194 | * MONIOU,MDEBUG,NUCNUC,
|
|---|
| 195 | * CETAPE,
|
|---|
| 196 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 197 | * N1STTR,MDBASE,
|
|---|
| 198 | * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 199 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 200 | * ,GHEISH,GHESIG
|
|---|
| 201 | COMMON /RUNPAC/ DSN,HOST,USER
|
|---|
| 202 | DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
|
|---|
| 203 | REAL STEPFC
|
|---|
| 204 | INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
|
|---|
| 205 | * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
|
|---|
| 206 | * N1STTR,MDBASE
|
|---|
| 207 | INTEGER CETAPE
|
|---|
| 208 | CHARACTER*79 DSN
|
|---|
| 209 | CHARACTER*20 HOST,USER
|
|---|
| 210 |
|
|---|
| 211 | LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
|
|---|
| 212 | * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
|
|---|
| 213 | * ,GHEISH,GHESIG
|
|---|
| 214 | *KEEP,STACKF.
|
|---|
| 215 | COMMON /STACKF/ STACK,STACKP,EXST,NSHIFT,NOUREC,ICOUNT,NTO,NFROM
|
|---|
| 216 | INTEGER MAXSTK
|
|---|
| 217 | PARAMETER (MAXSTK = 12*340*2)
|
|---|
| 218 | DOUBLE PRECISION STACK(MAXSTK)
|
|---|
| 219 | INTEGER STACKP,EXST,NSHIFT,NOUREC,ICOUNT,NTO,NFROM
|
|---|
| 220 | *KEEP,STATI.
|
|---|
| 221 | COMMON /STATI/ SABIN,SBBIN,INBIN,IPBIN,IKBIN,IHBIN
|
|---|
| 222 | DOUBLE PRECISION SABIN(37),SBBIN(37)
|
|---|
| 223 | INTEGER INBIN(37),IPBIN(37),IKBIN(37),IHBIN(37)
|
|---|
| 224 | *KEEP,THNVAR.
|
|---|
| 225 | COMMON /THNVAR/ STACKINT,INT_ICOUNT,THINNING
|
|---|
| 226 | INTEGER MAXICOUNT
|
|---|
| 227 | PARAMETER (MAXICOUNT=20000)
|
|---|
| 228 | DOUBLE PRECISION STACKINT(MAXICOUNT,13)
|
|---|
| 229 | INTEGER INT_ICOUNT
|
|---|
| 230 | LOGICAL THINNING
|
|---|
| 231 | *KEEP,VERS.
|
|---|
| 232 | COMMON /VERS/ VERNUM,MVDATE,VERDAT
|
|---|
| 233 | DOUBLE PRECISION VERNUM
|
|---|
| 234 | INTEGER MVDATE
|
|---|
| 235 | CHARACTER*18 VERDAT
|
|---|
| 236 | *KEEP,CEREN1.
|
|---|
| 237 | COMMON /CEREN1/ CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD,
|
|---|
| 238 | * CERSIZ,LCERFI
|
|---|
| 239 | DOUBLE PRECISION CERELE,CERHAD,ETADSN,WAVLGL,WAVLGU,CYIELD
|
|---|
| 240 | REAL CERSIZ
|
|---|
| 241 | LOGICAL LCERFI
|
|---|
| 242 | *KEEP,CEREN2.
|
|---|
| 243 | COMMON /CEREN2/ PHOTCM,XCER,YCER,UEMIS,VEMIS,CARTIM,ZEMIS,
|
|---|
| 244 | * DCERX,DCERY,ACERX,ACERY,
|
|---|
| 245 | * XCMAX,YCMAX,EPSX,EPSY,
|
|---|
| 246 | * DCERXI,DCERYI,FCERX,FCERY,
|
|---|
| 247 | * XSCATT,YSCATT,CERXOS,CERYOS,
|
|---|
| 248 | * NCERX,NCERY,ICERML
|
|---|
| 249 | REAL PHOTCM,XCER,YCER,UEMIS,VEMIS,CARTIM,ZEMIS,
|
|---|
| 250 | * DCERX,DCERY,ACERX,ACERY,
|
|---|
| 251 | * XCMAX,YCMAX,EPSX,EPSY,
|
|---|
| 252 | * DCERXI,DCERYI,FCERX,FCERY,
|
|---|
| 253 | * XSCATT,YSCATT,CERXOS(20),CERYOS(20)
|
|---|
| 254 | INTEGER NCERX,NCERY,ICERML
|
|---|
| 255 | *KEEP,CEREN3.
|
|---|
| 256 | COMMON /CEREN3/ CERCNT,DATAB2,LHCER
|
|---|
| 257 | INTEGER MAXBF2
|
|---|
| 258 | PARAMETER (MAXBF2 = 39 * 7)
|
|---|
| 259 | DOUBLE PRECISION CERCNT
|
|---|
| 260 | REAL DATAB2(MAXBF2)
|
|---|
| 261 | INTEGER LHCER
|
|---|
| 262 | *KEEP,CEREN4.
|
|---|
| 263 | COMMON /CEREN4/ NRECER
|
|---|
| 264 | INTEGER NRECER
|
|---|
| 265 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 266 | *keep,certel.
|
|---|
| 267 | common /certel/ cormxd,cord,coralp,ctpars,omega,
|
|---|
| 268 | + photn,photnp,phpt,pht,vphot,
|
|---|
| 269 | + vchi,veta,vzeta,vchim,vetam,vzetam,
|
|---|
| 270 | + lambda,mu,nu,nctels,ncph
|
|---|
| 271 | double precision cormxd,cord,coralp,ctpars(20,7),omega(20,3,3),
|
|---|
| 272 | + photn(3),photnp(3),phpt(3),pht,vphot(3),
|
|---|
| 273 | + vchi(3),veta(3),vzeta(3),vchim,vetam,vzetam,
|
|---|
| 274 | + lambda,mu,nu
|
|---|
| 275 | integer nctels,ncph(5)
|
|---|
| 276 | double precision xg,yg,zg,xgp,ygp,zgp,up,vp,wp,xpcut,ypcut,zpcut
|
|---|
| 277 | equivalence (photn(1) ,xg) ,(photn(2) ,yg) ,(photn(3) ,zg) ,
|
|---|
| 278 | + (photnp(1),xgp) ,(photnp(2),ygp) ,(photnp(3),zgp),
|
|---|
| 279 | + (phpt(1) ,xpcut),(phpt(2) ,ypcut),(phpt(3) ,zpcut),
|
|---|
| 280 | + (vphot(1) ,up) ,(vphot(2) ,vp) ,(vphot(3) ,wp)
|
|---|
| 281 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 282 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 283 | C Modificacion hecha por Aitor (5-feb-98)
|
|---|
| 284 | common /aitor/ aitoth
|
|---|
| 285 | double precision aitoth
|
|---|
| 286 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 287 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 288 | c Angles for the "spinning" of a particle around the
|
|---|
| 289 | c main axis of the CT
|
|---|
| 290 | common /spinang/ spinxi
|
|---|
| 291 | double precision spinxi
|
|---|
| 292 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 293 |
|
|---|
| 294 |
|
|---|
| 295 | *KEND.
|
|---|
| 296 |
|
|---|
| 297 | INTEGER JNBIN(37),JPBIN(37),JKBIN(37),JHBIN(37)
|
|---|
| 298 | DOUBLE PRECISION CHI2,FPARAM(6)
|
|---|
| 299 | DOUBLE PRECISION MPART2(10,25),MPHOT2(10),MPOSI2(10),MELEC2(10),
|
|---|
| 300 | * MNU2(10),MMUP2(10),MMUM2(10),MPI02(10),MPIP2(10),
|
|---|
| 301 | * MPIM2(10),MK0L2(10),MKPL2(10),MKMI2(10),
|
|---|
| 302 | * MNETR2(10),MPROT2(10),MPRTB2(10),MK0S2(10),
|
|---|
| 303 | * MHYP2(10),MNETB2(10),MDEUT2(10),MTRIT2(10),
|
|---|
| 304 | * MALPH2(10),MOTH2(10)
|
|---|
| 305 | EQUIVALENCE (MPART2(1, 1),MPHOT2(1)), (MPART2(1, 2),MPOSI2(1)),
|
|---|
| 306 | * (MPART2(1, 3),MELEC2(1)), (MPART2(1, 4),MNU2(1)) ,
|
|---|
| 307 | * (MPART2(1, 5),MMUP2(1)) , (MPART2(1, 6),MMUM2(1)) ,
|
|---|
| 308 | * (MPART2(1, 7),MPI02(1)) , (MPART2(1, 8),MPIP2(1)) ,
|
|---|
| 309 | * (MPART2(1, 9),MPIM2(1)) , (MPART2(1,10),MK0L2(1)) ,
|
|---|
| 310 | * (MPART2(1,11),MKPL2(1)) , (MPART2(1,12),MKMI2(1)) ,
|
|---|
| 311 | * (MPART2(1,13),MNETR2(1)), (MPART2(1,14),MPROT2(1)),
|
|---|
| 312 | * (MPART2(1,15),MPRTB2(1)), (MPART2(1,16),MK0S2(1)) ,
|
|---|
| 313 | * (MPART2(1,18),MHYP2(1)) , (MPART2(1,19),MDEUT2(1)),
|
|---|
| 314 | * (MPART2(1,20),MTRIT2(1)), (MPART2(1,21),MALPH2(1)),
|
|---|
| 315 | * (MPART2(1,22),MOTH2 (1)), (MPART2(1,25),MNETB2(1))
|
|---|
| 316 | C VARIABLES BEING USED FOR RUNTIME
|
|---|
| 317 | REAL TDIFF
|
|---|
| 318 | INTEGER ILEFTA,ILEFTB,TIME
|
|---|
| 319 | EXTERNAL TIME
|
|---|
| 320 |
|
|---|
| 321 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 322 | double precision ctdiams(20)
|
|---|
| 323 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 324 | double precision theprim, phiprim
|
|---|
| 325 | double precision spinthe, spinphi
|
|---|
| 326 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 327 |
|
|---|
| 328 |
|
|---|
| 329 | C-----------------------------------------------------------------------
|
|---|
| 330 |
|
|---|
| 331 |
|
|---|
| 332 | CERELE = 0.D0
|
|---|
| 333 | CERHAD = 0.D0
|
|---|
| 334 | NRECER = 0
|
|---|
| 335 | C INITIALIZE AND READ RUN STEERING CARDS
|
|---|
| 336 | CALL START
|
|---|
| 337 |
|
|---|
| 338 | IF ( CERSIZ .LE. 0. ) THEN
|
|---|
| 339 | ICRSIZ = 0
|
|---|
| 340 | ELSE
|
|---|
| 341 | ICRSIZ = 1
|
|---|
| 342 | ENDIF
|
|---|
| 343 |
|
|---|
| 344 | C RESET COUNTER FOR WORDS WRITTEN TO TAPE
|
|---|
| 345 | IRECOR = 0
|
|---|
| 346 |
|
|---|
| 347 | C RESET COUNTER FOR AVERAGE HIGHT OF 1ST INTERACTION
|
|---|
| 348 | CHISUM = 0.D0
|
|---|
| 349 | CHISM2 = 0.D0
|
|---|
| 350 |
|
|---|
| 351 | C SET ARRAYS FOR SCALES OF KINETIC ENERGY-INTERACTION TABLE
|
|---|
| 352 | SABIN(1) = 0.D0
|
|---|
| 353 | SBBIN(1) = 0.1D0
|
|---|
| 354 | DO 13 J = 2,37
|
|---|
| 355 | SABIN(J) = 10.D0**((J-4.D0)/3.D0)
|
|---|
| 356 | SBBIN(J) = 10.D0**((J-3.D0)/3.D0)
|
|---|
| 357 | 13 CONTINUE
|
|---|
| 358 |
|
|---|
| 359 | C CHECK AND SET PRIMARY PARAMETERS
|
|---|
| 360 | CALL INPRM
|
|---|
| 361 |
|
|---|
| 362 | do 161 i=1,nctels
|
|---|
| 363 | ctdiams(i) = ctpars(i,ctdiam)
|
|---|
| 364 | 161 continue
|
|---|
| 365 |
|
|---|
| 366 | C INITIALIZE NKG ROUTINES
|
|---|
| 367 | CALL ININKG
|
|---|
| 368 |
|
|---|
| 369 | C RESET COUNTERS FOR NUCLEON, PION AND KAON TABLE FOR ALL SHOWERS
|
|---|
| 370 | C RESET ENERGY-MULTIPLICITY & ENERGY-ELASTICITY MATRIX FOR ALL SHOWERS
|
|---|
| 371 | DO 17 J = 1,37
|
|---|
| 372 | JNBIN(J) = 0
|
|---|
| 373 | JPBIN(J) = 0
|
|---|
| 374 | JKBIN(J) = 0
|
|---|
| 375 | JHBIN(J) = 0
|
|---|
| 376 | ELMEAA(J) = 0.D0
|
|---|
| 377 | DO 17 L = 1,13
|
|---|
| 378 | MULTOT(J,L) = 0
|
|---|
| 379 | IELDPA(J,L) = 0
|
|---|
| 380 | 17 CONTINUE
|
|---|
| 381 |
|
|---|
| 382 | C RESET OTHER ARRAYS FOR ALL SHOWERS
|
|---|
| 383 | DO 99 J = 1,20
|
|---|
| 384 | IELNU(J) = 0
|
|---|
| 385 | IELPI(J) = 0
|
|---|
| 386 | IELIS(J) = 0
|
|---|
| 387 | IELHM(J) = 0
|
|---|
| 388 | 99 CONTINUE
|
|---|
| 389 |
|
|---|
| 390 | C RESET ARRAYS FOR INTERACTION LENGTH STATISTICS
|
|---|
| 391 | DO 90 J = 1,124
|
|---|
| 392 | IHYCHI(J) = 0
|
|---|
| 393 | IKACHI(J) = 0
|
|---|
| 394 | IMUCHI(J) = 0
|
|---|
| 395 | INUCHI(J) = 0
|
|---|
| 396 | IPICHI(J) = 0
|
|---|
| 397 | INNCHI(J) = 0
|
|---|
| 398 | 90 CONTINUE
|
|---|
| 399 |
|
|---|
| 400 | C RESET ARRAY FOR MEAN VALUES AND STANDARD DEVIATION
|
|---|
| 401 | DO 477 K = 1,25
|
|---|
| 402 | DO 477 J = 1,10
|
|---|
| 403 | MPARTO(J,K) = 0.D0
|
|---|
| 404 | MPART2(J,K) = 0.D0
|
|---|
| 405 | 477 CONTINUE
|
|---|
| 406 |
|
|---|
| 407 | C RESET ARRAYS FOR LONGITUDINAL DISTRIBUTION
|
|---|
| 408 | IF ( LLONGI ) THEN
|
|---|
| 409 | DO 478 K = 1,9
|
|---|
| 410 | DO 4781 J = 0,NSTEP
|
|---|
| 411 | APLONG(J,K) = 0.D0
|
|---|
| 412 | SPLONG(J,K) = 0.D0
|
|---|
| 413 | 4781 CONTINUE
|
|---|
| 414 | 478 CONTINUE
|
|---|
| 415 | ENDIF
|
|---|
| 416 |
|
|---|
| 417 | C STEERING OF PRINTOUT OF RANDOM GENERATOR SEEDS
|
|---|
| 418 | IPROUT = MIN(100,NSHOW/20)
|
|---|
| 419 | IPROUT = MAX(1,IPROUT)
|
|---|
| 420 |
|
|---|
| 421 | C TIME AT BEGINNING
|
|---|
| 422 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 423 | c ILEFTA = TIME()
|
|---|
| 424 | ILEFTA = 0
|
|---|
| 425 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 426 |
|
|---|
| 427 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 428 | print *,'JCIO::========================================'
|
|---|
| 429 | print *,'JCIO:: Initializing JCIO system for advanced'
|
|---|
| 430 | print *,'JCIO:: saving of data.'
|
|---|
| 431 | print *,'JCIO::========================================'
|
|---|
| 432 | c- initialize jcio system
|
|---|
| 433 | call jcinitio(dsn,nrrun)
|
|---|
| 434 | c- create file run######
|
|---|
| 435 | call jcstartrun(runh)
|
|---|
| 436 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 437 |
|
|---|
| 438 | C-----------------------------------------------------------------------
|
|---|
| 439 | C LOOP OVER SHOWERS
|
|---|
| 440 | DO 2 ISHW = 1,NSHOW
|
|---|
| 441 |
|
|---|
| 442 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 443 | c Next block of code has been modified, and comes from INPRM
|
|---|
| 444 | c----------------------------------------------------------------------
|
|---|
| 445 | C SCATTERING OF CENTER OF CHERENKOV ARRAY RELATIVE TO SHOWER AXIS
|
|---|
| 446 | c>> Actually, XSCATT and YSCATT should be RminSCATT and RmaxSCATT
|
|---|
| 447 | ICERML = MIN(MAX(ICERML,1),20)
|
|---|
| 448 | XSCATT = ABS(XSCATT)
|
|---|
| 449 | YSCATT = ABS(YSCATT)
|
|---|
| 450 | WRITE(MONIOU,5225)ICERML,XSCATT,YSCATT
|
|---|
| 451 | 5225 FORMAT(' ** USING EACH SHOWER SEVERAL TIMES:'/
|
|---|
| 452 | + ' USE EACH EVENT ',I2,' TIMES'/
|
|---|
| 453 | + ' THE EVENTS ARE SCATTERED RANDOMLY IN A SECTOR OF RADII:'/
|
|---|
| 454 | + ' Rmin = ',F10.0,' Rmax = ',F10.0)
|
|---|
| 455 | DO 4438 I=1,ICERML
|
|---|
| 456 | 5226 CALL RMMAR( RD,2,3 )
|
|---|
| 457 | CERXOS(I) = 2.0*YSCATT*(RD(1)-0.5)
|
|---|
| 458 | CERYOS(I) = 2.0*YSCATT*(RD(2)-0.5)
|
|---|
| 459 | R=SQRT(CERXOS(I)**2+CERYOS(I)**2)
|
|---|
| 460 | IF ((R.LT.XSCATT).OR.(R.GT.YSCATT)) GOTO 5226
|
|---|
| 461 | WRITE(MONIOU,4437) I,CERXOS(I),CERYOS(I)
|
|---|
| 462 | 4437 FORMAT(' CORE OF EVENT ',I2,' AT ',2F12.2)
|
|---|
| 463 | 4438 CONTINUE
|
|---|
| 464 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 465 |
|
|---|
| 466 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 467 | c Next block of code comes from INPRM
|
|---|
| 468 | c----------------------------------------------------------------------
|
|---|
| 469 | EVTH(98) = FLOAT(ICERML)
|
|---|
| 470 | DO 480 I=1,20
|
|---|
| 471 | EVTH( 98+I) = CERXOS(I)
|
|---|
| 472 | EVTH(118+I) = CERYOS(I)
|
|---|
| 473 | 480 CONTINUE
|
|---|
| 474 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 475 |
|
|---|
| 476 | SHOWNO = SHOWNO + 1
|
|---|
| 477 | I = ISHW
|
|---|
| 478 | IF ( ISHW .LE. MAXPRT ) THEN
|
|---|
| 479 | FPRINT = .TRUE.
|
|---|
| 480 | ELSE
|
|---|
| 481 | FPRINT = .FALSE.
|
|---|
| 482 | ENDIF
|
|---|
| 483 |
|
|---|
| 484 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 485 | c Create cer######,dat######,sta###### files
|
|---|
| 486 | c------------------------------------------------------------
|
|---|
| 487 | call jcnewshower
|
|---|
| 488 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 489 |
|
|---|
| 490 | C RESET COUNTERS
|
|---|
| 491 | DO 447 K = 1,25
|
|---|
| 492 | DO 447 J = 1,10
|
|---|
| 493 | NPARTO(J,K) = 0.D0
|
|---|
| 494 | 447 CONTINUE
|
|---|
| 495 | MUOND = 0.D0
|
|---|
| 496 |
|
|---|
| 497 | C RESET ARRAY FOR LONGITUDINAL DISTRIBUTION PER SHOWER
|
|---|
| 498 | IF ( LLONGI ) THEN
|
|---|
| 499 | DO 479 K = 1,9
|
|---|
| 500 | DO 4791 J = 0,NSTEP
|
|---|
| 501 | PLONG(J,K) = 0.D0
|
|---|
| 502 | 4791 CONTINUE
|
|---|
| 503 | 479 CONTINUE
|
|---|
| 504 | ENDIF
|
|---|
| 505 |
|
|---|
| 506 | NRECS = 0
|
|---|
| 507 | NBLKS = 0
|
|---|
| 508 | DO 922 KKK = 1,10
|
|---|
| 509 | AVNREJ(KKK) = 0.D0
|
|---|
| 510 | 922 CONTINUE
|
|---|
| 511 | IRESPAR = 0
|
|---|
| 512 |
|
|---|
| 513 |
|
|---|
| 514 | C FIRST INTERACTION DATA
|
|---|
| 515 | FIRSTI = .TRUE.
|
|---|
| 516 | IFINET = 0
|
|---|
| 517 | IFINNU = 0
|
|---|
| 518 | IFINKA = 0
|
|---|
| 519 | IFINPI = 0
|
|---|
| 520 | IFINHY = 0
|
|---|
| 521 | ELAST = 0.D0
|
|---|
| 522 |
|
|---|
| 523 | C RESET COUNTERS FOR NUCLEON, PION AND KAON TABLE FOR SHOWER
|
|---|
| 524 | C RESET ENERGY-MULTIPLICITY & ENERGY-ELASTICITY MATRIX FOR SHOWER
|
|---|
| 525 | DO 11 J = 1,37
|
|---|
| 526 | INBIN(J) = 0
|
|---|
| 527 | IPBIN(J) = 0
|
|---|
| 528 | IKBIN(J) = 0
|
|---|
| 529 | IHBIN(J) = 0
|
|---|
| 530 | ELMEAN(J) = 0.D0
|
|---|
| 531 | DO 11 L = 1,13
|
|---|
| 532 | MULTMA(J,L) = 0
|
|---|
| 533 | IELDPM(J,L) = 0
|
|---|
| 534 | 11 CONTINUE
|
|---|
| 535 |
|
|---|
| 536 | DO 12 J = 1,10
|
|---|
| 537 | PBAL(J) = 0.D0
|
|---|
| 538 | EBAL(J) = 0.D0
|
|---|
| 539 | 12 CONTINUE
|
|---|
| 540 |
|
|---|
| 541 | C INITIALIZE PARTICLE STACK
|
|---|
| 542 | CALL ISTACK
|
|---|
| 543 | C RESET STACKINT
|
|---|
| 544 | DO J=1,MAXICOUNT
|
|---|
| 545 | DO K=1,MAXLEN
|
|---|
| 546 | STACKINT(J,K) = 0.D0
|
|---|
| 547 | ENDDO
|
|---|
| 548 | ENDDO
|
|---|
| 549 |
|
|---|
| 550 | C INITIALIZE EVENT HEADER AND END FOR EACH EVENT
|
|---|
| 551 | DO 2123 L = 2,43
|
|---|
| 552 | EVTH(L) = 0.
|
|---|
| 553 | 2123 CONTINUE
|
|---|
| 554 | DO 123 L = 2,MAXBUF
|
|---|
| 555 | EVTE(L) = 0.
|
|---|
| 556 | 123 CONTINUE
|
|---|
| 557 |
|
|---|
| 558 | C SHOWER BEGIN PRINTOUT
|
|---|
| 559 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,105) SHOWNO
|
|---|
| 560 | 105 FORMAT ('1',10('='),' SHOWER NO ',I10,' ',47('=')/)
|
|---|
| 561 |
|
|---|
| 562 | C RANDOM GENERATOR STATUS AT BEGINNING OF SHOWER CALCULATION
|
|---|
| 563 | EVTH(13) = NSEQ
|
|---|
| 564 | DO 45 L = 1,NSEQ
|
|---|
| 565 | CALL RMMAQ( ISEED(1,L), L, 'R' )
|
|---|
| 566 | C SEED
|
|---|
| 567 | EVTH(11+L*3) = ISEED(1,L)
|
|---|
| 568 | C NUMBER OF CALLS
|
|---|
| 569 | EVTH(12+L*3) = MOD ( ISEED(2,L), 1000000 )
|
|---|
| 570 | C NUMBER OF MILLIONS
|
|---|
| 571 | EVTH(13+L*3) = ISEED(3,L)*1000 + INT( ISEED(2,L)/1000000 )
|
|---|
| 572 | 45 CONTINUE
|
|---|
| 573 | IF ( FPRINT .OR. DEBUG .OR. MOD(ISHW-1,IPROUT).EQ.0 ) THEN
|
|---|
| 574 | CALL PRTIME(TTIME)
|
|---|
| 575 | WRITE(MONIOU,158) SHOWNO,(L,(ISEED(J,L),J=1,3),L=1,NSEQ)
|
|---|
| 576 | 158 FORMAT(' AND RANDOM NUMBER GENERATOR AT BEGIN OF EVENT :',I8,
|
|---|
| 577 | * /,(' SEQUENCE = ',I2,' SEED = ',I9 ,' CALLS = ',I9,
|
|---|
| 578 | * ' BILLIONS = ',I9))
|
|---|
| 579 | ENDIF
|
|---|
| 580 | C RESET KNOR
|
|---|
| 581 | KNOR = .TRUE.
|
|---|
| 582 |
|
|---|
| 583 | C GET FULL RANDOM GENERATOR STATUS (103 WORDS PER SEQUENCE)
|
|---|
| 584 | CC DO 495 L = 1,NSEQ
|
|---|
| 585 | CC CALL RMMAQ( ISEED(1,L), L, 'RV' )
|
|---|
| 586 | CC WRITE(MONIOU,658) L,(ISEED(J,L),J=1,103)
|
|---|
| 587 | CC658 FORMAT ( ' FULL RANDOM NUMBER GENERATOR STATUS ',
|
|---|
| 588 | CC * 'FOR SEQUENCE ',I2,/(' ',10I11))
|
|---|
| 589 | CC495 CONTINUE
|
|---|
| 590 |
|
|---|
| 591 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 592 | c>> *** ATENTION *** ATENTION *** ATENTION *** ATENTION *** ATENTION >>
|
|---|
| 593 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 594 | c>> >>
|
|---|
| 595 | c>> In the next block (between this ATENTION comments) CORSIKA makes >>
|
|---|
| 596 | c>> three things, in this order: >>
|
|---|
| 597 | c>> >>
|
|---|
| 598 | c>> i. Set ANGLES OF INCIDENCE (different distributions of theta >>
|
|---|
| 599 | c>> for gammas -flat- and hadrons -standard. >>
|
|---|
| 600 | c>> ii. Set HEIGHT for start at THICK0 (normally = 0 => 112.8 Km) >>
|
|---|
| 601 | c>> iii. Set ENERGY of the primary. >>
|
|---|
| 602 | c>> >>
|
|---|
| 603 | c>> (The original order was ii., iii. and i.) >>
|
|---|
| 604 | c>> >>
|
|---|
| 605 | c>> JCG Wed Sep 21 10:49:14 MET DST 1998 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 606 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 607 |
|
|---|
| 608 |
|
|---|
| 609 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 610 | C GET PRIMARY ANGLES OF INCIDENCE
|
|---|
| 611 | IF ( FIXINC ) THEN
|
|---|
| 612 |
|
|---|
| 613 | THETAP = THETPR(1)
|
|---|
| 614 | PHIP = PHIPR(1)
|
|---|
| 615 | PRMPAR(3) = COS(THETAP)
|
|---|
| 616 |
|
|---|
| 617 | ELSE
|
|---|
| 618 |
|
|---|
| 619 | if ( prmpar(1).eq.1 ) then
|
|---|
| 620 |
|
|---|
| 621 | C>> GAMMAS >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 622 | C NOTE!! We will use a FLAT distribution for THETA:
|
|---|
| 623 | C Then, next block (original block) must be commented.
|
|---|
| 624 | c The modificated code follows this block
|
|---|
| 625 | c
|
|---|
| 626 | CALL RMMAR( RD,3,1 )
|
|---|
| 627 | CT1 = THETPR(1)
|
|---|
| 628 | CT2 = THETPR(2)
|
|---|
| 629 | THETAP = RD(2)*(CT2 - CT1) + CT1
|
|---|
| 630 | CTT = COS(THETAP)
|
|---|
| 631 | PRMPAR(3) = CTT
|
|---|
| 632 |
|
|---|
| 633 | else
|
|---|
| 634 |
|
|---|
| 635 | C>> HADRONS AND ELECTRONS (AND ANY OTHER BUT GAMMAS) >>>>>>>>>>>>>>>>
|
|---|
| 636 | c Choose angles at random with equal flux for all directions
|
|---|
| 637 | c with horizontal detector array (see: O.C. Allkofer & P.K.F. Grieder,
|
|---|
| 638 | c Cosmic Rays on Earth, in: Physics Data 25/1, H.Behrens & G.Ebel Ed.,
|
|---|
| 639 | c (Fachinformationszentrum Karlsruhe, Germany, 1983) chpt. 1.1.2)
|
|---|
| 640 | c
|
|---|
| 641 | CALL RMMAR( RD,3,1 )
|
|---|
| 642 | CT1 = SIN(THETPR(1))**2
|
|---|
| 643 | CT2 = SIN(THETPR(2))**2
|
|---|
| 644 | CTT = SQRT( 1.D0 - RD(2)*(CT2 - CT1) - CT1 )
|
|---|
| 645 | PRMPAR(3) = CTT
|
|---|
| 646 | THETAP = ACOS(CTT)
|
|---|
| 647 |
|
|---|
| 648 | endif
|
|---|
| 649 |
|
|---|
| 650 | PHIP = RD(1) * ( PHIPR(2) - PHIPR(1) ) + PHIPR(1)
|
|---|
| 651 |
|
|---|
| 652 | ENDIF
|
|---|
| 653 | PRMPAR(4) = PHIP
|
|---|
| 654 |
|
|---|
| 655 |
|
|---|
| 656 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 657 | C DEFINE HEIGHT FOR START AT THICK0 (IN G/CM**2) (112.8 KM FOR THICK0=0)
|
|---|
| 658 | PRMPAR(5) = HEIGH(THICK0)
|
|---|
| 659 |
|
|---|
| 660 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 661 | C GET PRIMARY ENERGY INTO PRMPAR(2)
|
|---|
| 662 | IF ( ISPEC .EQ. 0 ) THEN
|
|---|
| 663 | PRMPAR(2) = LLIMIT
|
|---|
| 664 | ELSE
|
|---|
| 665 | CALL RMMAR( RD,1,1 )
|
|---|
| 666 | IF ( PSLOPE .NE. -1.D0 ) THEN
|
|---|
| 667 | PRMPAR(2) = ( RD(1)*UL + ( 1.D0-RD(1) )*LL )**SLEX
|
|---|
| 668 | ELSE
|
|---|
| 669 | PRMPAR(2) = LLIMIT * LL**RD(1)
|
|---|
| 670 | ENDIF
|
|---|
| 671 | ENDIF
|
|---|
| 672 |
|
|---|
| 673 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 674 | c>> *** ATENTION *** ATENTION *** ATENTION *** ATENTION *** ATENTION >>
|
|---|
| 675 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 676 |
|
|---|
| 677 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 678 | c>> Modification: this is no longer needed >>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 679 | c>> (Superseeded by Sphere algorithm, see cerenkov.f) >>>>>>>>>>>>>>>>>
|
|---|
| 680 | c>> JCG Wed Sep 21 10:49:14 MET DST 1998 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 681 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 682 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 683 | cc Btw, we now modify the "shadow area" of the telescopes,
|
|---|
| 684 | cc to cover the angle theta.
|
|---|
| 685 | c do 160 i=1,nctels
|
|---|
| 686 | c ctpars(i,ctdiam) = ctdiams(i)/cos(thetap)
|
|---|
| 687 | c write (MONIOU,*)
|
|---|
| 688 | c * 'New region for CT',i,' = ',ctpars(i,ctdiam)
|
|---|
| 689 | c 160 continue
|
|---|
| 690 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 691 |
|
|---|
| 692 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 693 | c Here we calculate the angles spinphi and spinthe, which are the
|
|---|
| 694 | c phi and theta angles of the particle, with respect to the direction
|
|---|
| 695 | c where the CT is pointing to. spinthe is the angular displacement
|
|---|
| 696 | c of the new direction respect to the original (CT); spinphi=0 means
|
|---|
| 697 | c that the new direction is towards the zenith, spinphi=+-180 means
|
|---|
| 698 | c towards the horizont.
|
|---|
| 699 | c See the document "simulation.tex"
|
|---|
| 700 |
|
|---|
| 701 | c First, save the "CT" orientation
|
|---|
| 702 | c (moved from a couple of lines below, marked with [*])
|
|---|
| 703 | EVTH(11) = THETAP
|
|---|
| 704 | EVTH(12) = PHIP
|
|---|
| 705 |
|
|---|
| 706 | CALL RMMAR( RD,3,1 )
|
|---|
| 707 |
|
|---|
| 708 | c Then, calculate the new direction relative to the CT direction
|
|---|
| 709 | spinphi = RD(1)*PI
|
|---|
| 710 | spinthe = RD(2)*spinxi*pi/180
|
|---|
| 711 |
|
|---|
| 712 | c And then, RE-calculate the GLOBAL direction in CORSIKA
|
|---|
| 713 | c We use formulae for spherical triangles
|
|---|
| 714 | theprim = acos( cos(THETAP)*cos(spinthe)+
|
|---|
| 715 | $ sin(THETAP)*sin(spinthe)*cos(spinphi) )
|
|---|
| 716 | phiprim = asin( sin(spinthe)*sin(spinphi)/sin(theprim) )
|
|---|
| 717 | THETAP = theprim
|
|---|
| 718 | EVTH(140) = spinthe
|
|---|
| 719 | if (RD(3).gt.0.5) then
|
|---|
| 720 | PHIP = PHIP - phiprim
|
|---|
| 721 | EVTH(141) = -spinphi
|
|---|
| 722 | else
|
|---|
| 723 | PHIP = PHIP + phiprim
|
|---|
| 724 | EVTH(141) = spinphi
|
|---|
| 725 | endif
|
|---|
| 726 | PRMPAR(3) = COS(THETAP)
|
|---|
| 727 | PRMPAR(4) = PHIP
|
|---|
| 728 |
|
|---|
| 729 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 730 |
|
|---|
| 731 |
|
|---|
| 732 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 733 | c>> Modification (HZA trick) cancelled >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 734 | c>> JCG Wed Sep 21 10:49:14 MET DST 1998 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 735 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 736 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 737 | C Modificacion hecha por Aitor
|
|---|
| 738 | c aitoth = THETAP
|
|---|
| 739 | C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 740 |
|
|---|
| 741 | C CALCULATE COORDINATE CORRECTION FOR TOP OF ATMOSPHERE
|
|---|
| 742 | CALL COORIN( PRMPAR(5) )
|
|---|
| 743 |
|
|---|
| 744 | C COUNTER FOR PARTICLE OUTPUT
|
|---|
| 745 | LH = 0
|
|---|
| 746 | C COUNTER FOR CERENKOV OUTPUT
|
|---|
| 747 | IF ( LCERFI ) LHCER = 0
|
|---|
| 748 | C CALCULATE BUNCH SIZE FOR CERENKOV PHOTONS IF NOT SET IN DATAC
|
|---|
| 749 | IF ( ICRSIZ .EQ. 0 ) THEN
|
|---|
| 750 | CALL GETBUS( NINT(PRMPAR(1)),SNGL(PRMPAR(2)),SNGL(PRMPAR(3)),
|
|---|
| 751 | * CERSIZ )
|
|---|
| 752 | WRITE(MONIOU,*)'CERENKOV BUNCH SIZE IS CALCULATED TO=',CERSIZ
|
|---|
| 753 | ENDIF
|
|---|
| 754 | C GET GAMMA FACTOR FROM ENERGY
|
|---|
| 755 | C FOR GAMMAS PRMPAR(2) STAYS = ENERGY
|
|---|
| 756 | IF ( PRMPAR(1) .NE. 1.D0 )
|
|---|
| 757 | * PRMPAR(2) = PRMPAR(2) / PAMA(NINT(PRMPAR(1)))
|
|---|
| 758 |
|
|---|
| 759 | C SET PRIMARY TO CURRENT PARTICLE
|
|---|
| 760 | DO 3 J = 1,8
|
|---|
| 761 | CURPAR(J) = PRMPAR(J)
|
|---|
| 762 | NCOUN(J) = 0
|
|---|
| 763 | 3 CONTINUE
|
|---|
| 764 | C SET WEIGHT
|
|---|
| 765 |
|
|---|
| 766 | C CALCULATE FIRST INTERACTION POINT IF HADRONIC
|
|---|
| 767 | GEN = 0.D0
|
|---|
| 768 |
|
|---|
| 769 | H = HEIGH(THICK0)
|
|---|
| 770 | CALL BOX2
|
|---|
| 771 | IF ( FIX1I ) THEN
|
|---|
| 772 | CHI = THICK(FIXHEI) / PRMPAR(3)
|
|---|
| 773 | H = FIXHEI
|
|---|
| 774 | FDECAY = .FALSE.
|
|---|
| 775 | ELSE
|
|---|
| 776 | H = HEIGH ( CHI*PRMPAR(3) + THICK0 )
|
|---|
| 777 | ENDIF
|
|---|
| 778 | CHISUM = CHISUM + CHI
|
|---|
| 779 | CHISM2 = CHISM2 + CHI**2
|
|---|
| 780 | ALEVEL = H
|
|---|
| 781 | C INITIALIZE COORDINATE CORRECTIONS FOR HADRONIC PRIMARIES
|
|---|
| 782 | C FOR EM PRIMARIES IT IS DONE IN EGS
|
|---|
| 783 | HH = MAX( H, 0.D0 )
|
|---|
| 784 | IF ( CURPAR(1) .GT. 3.D0 ) CALL COORIN( HH )
|
|---|
| 785 |
|
|---|
| 786 | IF ( FMUADD ) THEN
|
|---|
| 787 | IF ( CURPAR(1) .EQ. 5 .OR. CURPAR(1) .EQ. 6) THEN
|
|---|
| 788 | DO J = 1,MAXLEN
|
|---|
| 789 | AMUPAR(J) = CURPAR(J)
|
|---|
| 790 | ENDDO
|
|---|
| 791 | AMUPAR(5) = PRMPAR(5)
|
|---|
| 792 | IF(DEBUG)WRITE(MDEBUG,*)'MAIN : MUON STORED IN AMUPAR'
|
|---|
| 793 | FMUORG = .TRUE.
|
|---|
| 794 | ENDIF
|
|---|
| 795 | ENDIF
|
|---|
| 796 |
|
|---|
| 797 | C SET TARGET FLAG IF SELECTED FOR FIRST INTERACTION
|
|---|
| 798 | IF ( N1STTR .GT. 0 ) THEN
|
|---|
| 799 | FIXTAR = .TRUE.
|
|---|
| 800 | FDECAY = .FALSE.
|
|---|
| 801 | EVTH(6) = REAL(N1STTR)
|
|---|
| 802 | ELSE
|
|---|
| 803 | FIXTAR = .FALSE.
|
|---|
| 804 | EVTH(6) = 0.0
|
|---|
| 805 | ENDIF
|
|---|
| 806 |
|
|---|
| 807 | C INITIALIZE ARRAYS FOR NKG FOR EACH SHOWER
|
|---|
| 808 | IF ( FNKG ) CALL STANKG
|
|---|
| 809 |
|
|---|
| 810 | C STORE FIRST PARTICLE IN HEADER AND PRINT IT OUT
|
|---|
| 811 | EVTH( 2) = REAL(SHOWNO)
|
|---|
| 812 | EVTH( 3) = CURPAR(1)
|
|---|
| 813 | IF ( CURPAR(1) .EQ. 1.D0 ) THEN
|
|---|
| 814 | C PRIMARY ENERGY FOR PHOTONS
|
|---|
| 815 | E00 = GAMMA
|
|---|
| 816 | E00PN = GAMMA
|
|---|
| 817 | INUCL = 1
|
|---|
| 818 | ELSE
|
|---|
| 819 | E00 = GAMMA * PAMA(NINT(CURPAR(1)))
|
|---|
| 820 | INUCL = INT(MAX(1.D0,CURPAR(1)/100.D0))
|
|---|
| 821 | E00PN = E00 / INUCL
|
|---|
| 822 | ENDIF
|
|---|
| 823 | EVTH(147) = 0.
|
|---|
| 824 |
|
|---|
| 825 | IF ( FEGS ) THEN
|
|---|
| 826 | C PARAMETER FOR ELECTRON AND PHOTON REJECT (CONVERT ENERGY TO MEV)
|
|---|
| 827 | EONCUT = .5E-9*SQRT(E00*1000.D0)
|
|---|
| 828 | CUTLN = LOG(EONCUT)
|
|---|
| 829 | ENDIF
|
|---|
| 830 | EVTH( 4) = E00
|
|---|
| 831 | EVTH( 5) = THICK0
|
|---|
| 832 | EVTH( 7) = H
|
|---|
| 833 | PTOT0 = SQRT( E00**2 - PAMA(NINT(CURPAR(1)))**2 )
|
|---|
| 834 | PTOT0N = PTOT0 / INUCL
|
|---|
| 835 | ST = SQRT(1.D0-COSTHE**2)
|
|---|
| 836 | EVTH( 8) = PTOT0 * ST * COS(PHI)
|
|---|
| 837 | EVTH( 9) = PTOT0 * ST * SIN(PHI)
|
|---|
| 838 | EVTH(10) = PTOT0 * COSTHE
|
|---|
| 839 | c
|
|---|
| 840 | c [*] one block from here sent above
|
|---|
| 841 | c
|
|---|
| 842 | EVTH(85) = CERSIZ
|
|---|
| 843 |
|
|---|
| 844 | IF ( CURPAR(1) .GT. 3.D0 ) THEN
|
|---|
| 845 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,102) (CURPAR(J),J = 1,8)
|
|---|
| 846 | 102 FORMAT (/' PRIMARY PARAMETERS AT FIRST INTERACTION POINT'/
|
|---|
| 847 | * 16X,1P,8E10.3)
|
|---|
| 848 | ELSE
|
|---|
| 849 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,132)
|
|---|
| 850 | 132 FORMAT (/' PRIMARY PARTICLE IS ELECTROMAGNETIC')
|
|---|
| 851 | ENDIF
|
|---|
| 852 |
|
|---|
| 853 | C WRITE EVENT HEADER INTO BUFFER
|
|---|
| 854 | C FOR EM PARTICLES EVTH IS WRITTEN TO BUFFER IN EGS (IF ACTIVE)
|
|---|
| 855 | IF ( EVTH(3) .GT. 3.0 .OR. .NOT. FEGS ) THEN
|
|---|
| 856 | CALL TOBUF ( EVTH,0 )
|
|---|
| 857 | IF ( LCERFI ) CALL TOBUFC( EVTH,0 )
|
|---|
| 858 | ENDIF
|
|---|
| 859 |
|
|---|
| 860 | C PRINT HEADER FOR HIGH ENERGY PARTICLES
|
|---|
| 861 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,103)
|
|---|
| 862 | 103 FORMAT(/' TYPE GAMMA COSTHETA ',
|
|---|
| 863 | * ' PHI HEIGHT TIME X-CM Y-CM ',
|
|---|
| 864 | * ' GEN LEVEL E ON STACK'/)
|
|---|
| 865 | NOPART = 0
|
|---|
| 866 |
|
|---|
| 867 |
|
|---|
| 868 | IF ( CURPAR(1) .LE. 3.D0 .OR.
|
|---|
| 869 | * (CURPAR(1) .GE. 5.D0 .AND. CURPAR(1) .LE. 7.D0) ) THEN
|
|---|
| 870 | C GIVE PARTICLE TO EGS OR NKG IF ELECTROMAGNETIC
|
|---|
| 871 | C AND TAKE THEN NEXT PARTICLE FROM STACK
|
|---|
| 872 | C FLAG FOR NO PRIMARY INTERACTION IS SET FOR ALL BUT ELM. PRIMARIES
|
|---|
| 873 | IF ( CURPAR(1) .LE. 3.D0 ) THEN
|
|---|
| 874 | FNPRIM = .FALSE.
|
|---|
| 875 | ELSE
|
|---|
| 876 | FNPRIM = .TRUE.
|
|---|
| 877 | H = PRMPAR(5)
|
|---|
| 878 | ENDIF
|
|---|
| 879 | CALL BOX3
|
|---|
| 880 | IF ( FEGS ) THEN
|
|---|
| 881 | CHISUM = CHISUM + THICK(DBLE(EVTH(7)))
|
|---|
| 882 | CHISM2 = CHISM2 + THICK(DBLE(EVTH(7)))**2
|
|---|
| 883 | ENDIF
|
|---|
| 884 | FIRSTI = .FALSE.
|
|---|
| 885 | GOTO 4
|
|---|
| 886 |
|
|---|
| 887 | ELSE
|
|---|
| 888 | C HADRONIC PARTICLES
|
|---|
| 889 | FNPRIM = .TRUE.
|
|---|
| 890 |
|
|---|
| 891 | C FILL LONGITUDINAL DISTRIBUTION FOR THE PRIMARY PARTICLE
|
|---|
| 892 | C THE PARTICLE IS TRACKED FROM THICK0 DOWN TO THICK0+CHI*PRMPAR(3)
|
|---|
| 893 | C COUNT THE PARTICLES FOR THE LONGITUDINAL DEVELOPMENT
|
|---|
| 894 | IF ( LLONGI ) THEN
|
|---|
| 895 | LPCT1 = INT( THICK0 * THSTPI )
|
|---|
| 896 | LPCT2 = INT( (THICK0 + PRMPAR(3)*CHI) * THSTPI )
|
|---|
| 897 | LPCT2 = MIN(NSTEP,LPCT2)
|
|---|
| 898 | C GAMMAS, ELECTRONS AND POSITRONS ARE NOT TRANSPORTED HERE, SEE EGS
|
|---|
| 899 | C MUONS ARE TRANSPORTED IN MUTRAC
|
|---|
| 900 | C HADRONS
|
|---|
| 901 | IF ( ITYPE .GE. 7 .AND. ITYPE .LE. 41 ) THEN
|
|---|
| 902 | DO 5004 L = LPCT1,LPCT2
|
|---|
| 903 | PLONG(L,6) = PLONG(L,6) + 1.
|
|---|
| 904 | 5004 CONTINUE
|
|---|
| 905 | C CHARGED HADRONS
|
|---|
| 906 | IF ( SIGNUM(ITYPE) .NE. 0.D0 ) THEN
|
|---|
| 907 | DO 5005 L = LPCT1,LPCT2
|
|---|
| 908 | PLONG(L,7) = PLONG(L,7) + 1.
|
|---|
| 909 | 5005 CONTINUE
|
|---|
| 910 | ENDIF
|
|---|
| 911 | C NUCLEI
|
|---|
| 912 | ELSEIF ( ITYPE .GT. 100 ) THEN
|
|---|
| 913 | DO 5006 L = LPCT1,LPCT2
|
|---|
| 914 | PLONG(L,8) = PLONG(L,8) + 1.
|
|---|
| 915 | 5006 CONTINUE
|
|---|
| 916 | ENDIF
|
|---|
| 917 | ENDIF
|
|---|
| 918 |
|
|---|
| 919 | C CHECK OBSERVATION LEVEL PASSAGE AND UPDATE PARTICLE COORDINATES
|
|---|
| 920 | HNEW = H
|
|---|
| 921 | C FOR UPDATE WE NEED THE START ALTITUDE H
|
|---|
| 922 | H = HEIGH(THICK0)
|
|---|
| 923 | DO 251 J = 1,NOBSLV
|
|---|
| 924 | C JUMP INTO NORMAL PARTICLE TREATMENT FOR HADRONS
|
|---|
| 925 | IF ( HNEW .GT. OBSLEV(J) ) THEN
|
|---|
| 926 | H = HNEW
|
|---|
| 927 | GOTO 6
|
|---|
| 928 | ENDIF
|
|---|
| 929 | IF ( H .LT. OBSLEV(J) ) GOTO 251
|
|---|
| 930 | C REMEMBER NUMBER OF LEVEL FOR OUTPUT
|
|---|
| 931 | LEVL = J
|
|---|
| 932 | CALL UPDATE( OBSLEV(J),THCKOB(J),J )
|
|---|
| 933 | IF (DEBUG) WRITE(MDEBUG,256) J,IRET1,IRET2
|
|---|
| 934 | 256 FORMAT(' MAIN : LEVEL ',I5,' IRET1,2=',2I5)
|
|---|
| 935 | C IF PARTICLE IS NOT CUTTED, BRING IT TO OUTPUT
|
|---|
| 936 | IF ( IRET2 .EQ. 0 ) THEN
|
|---|
| 937 | CALL OUTPUT
|
|---|
| 938 | ENDIF
|
|---|
| 939 | 251 CONTINUE
|
|---|
| 940 | IF (DEBUG) WRITE(MDEBUG,*)
|
|---|
| 941 | * 'MAIN : PRIMARY REACHED LOWEST OBSERVATION LEVEL'
|
|---|
| 942 | GOTO 4
|
|---|
| 943 | ENDIF
|
|---|
| 944 |
|
|---|
| 945 | C-----------------------------------------------------------------------
|
|---|
| 946 | C NORMAL CYCLE
|
|---|
| 947 | 7 CONTINUE
|
|---|
| 948 |
|
|---|
| 949 | C IF ENERGY TOO SMALL TAKE NEXT PARTICLE
|
|---|
| 950 | IF ( GAMMA .LE. 1.D0 ) THEN
|
|---|
| 951 | IF ( CURPAR(1) .NE. 1.D0 ) THEN
|
|---|
| 952 | IF ( CURPAR(1).EQ.5.D0 .OR. CURPAR(1).EQ.6.D0 )
|
|---|
| 953 | * FMUORG = .FALSE.
|
|---|
| 954 | GOTO 4
|
|---|
| 955 | ENDIF
|
|---|
| 956 | C SPECIAL TREATMENT FOR PHOTONS
|
|---|
| 957 | ITYPE = 1
|
|---|
| 958 | CHI = 0.D0
|
|---|
| 959 | GOTO 5
|
|---|
| 960 | ENDIF
|
|---|
| 961 |
|
|---|
| 962 | C DETERMINE PLACE OF NEXT INTERACTION
|
|---|
| 963 | CALL BOX2
|
|---|
| 964 |
|
|---|
| 965 | C CHECK PASSAGE THROUGH OBSERVATION LEVELS AND TRACK PARTICLES TO THE
|
|---|
| 966 | C PLACE OF INTERACTION
|
|---|
| 967 | 5 CONTINUE
|
|---|
| 968 | IRET1 = 0
|
|---|
| 969 | CALL BOX3
|
|---|
| 970 | IF ( IRET1 .NE. 0 ) GOTO 4
|
|---|
| 971 |
|
|---|
| 972 | 6 CONTINUE
|
|---|
| 973 | IRET1 = 0
|
|---|
| 974 | MSMM = 0
|
|---|
| 975 |
|
|---|
| 976 | C INCREMENT PARTICLE GENERATION AND PROCESS NUCLEAR INTERACTION
|
|---|
| 977 | GEN = GEN + 1.D0
|
|---|
| 978 | C INITIALIZE INTERMEDIATE STACK FOR ONE REACTION
|
|---|
| 979 | CALL TSTINI
|
|---|
| 980 | CALL NUCINT
|
|---|
| 981 | C TRANSFER INTERMEDIATE STACK FOR ONE REACTION
|
|---|
| 982 | CALL TSTEND
|
|---|
| 983 |
|
|---|
| 984 | C ENERGY - MULTIPLICITY STATISTICS
|
|---|
| 985 | IF ( EKINL .LE. 0.1D0 ) THEN
|
|---|
| 986 | MEN = 1
|
|---|
| 987 | ELSE
|
|---|
| 988 | MEN = 4.D0 + 3.D0 * LOG10(EKINL)
|
|---|
| 989 | MEN = MIN( MEN, 37 )
|
|---|
| 990 | ENDIF
|
|---|
| 991 | IF ( MSMM .LE. 1 ) THEN
|
|---|
| 992 | MMU = 1
|
|---|
| 993 | ELSE
|
|---|
| 994 | MMU = 1.D0 + 3.D0 * LOG10(DBLE(MSMM))
|
|---|
| 995 | MMU = MIN( MMU, 13 )
|
|---|
| 996 | ENDIF
|
|---|
| 997 | MULTMA(MEN,MMU) = MULTMA(MEN,MMU) + 1
|
|---|
| 998 | MULTOT(MEN,MMU) = MULTOT(MEN,MMU) + 1
|
|---|
| 999 | IF ( DEBUG ) WRITE(MDEBUG,*) 'MAIN : EKINL,MSMM=',
|
|---|
| 1000 | * SNGL(EKINL),MSMM
|
|---|
| 1001 |
|
|---|
| 1002 | IF ( IRET1 .EQ. 0 ) THEN
|
|---|
| 1003 | IF ( DEBUG ) WRITE(MDEBUG,666) (CURPAR(II),II=1,11)
|
|---|
| 1004 | 666 FORMAT(' MAIN : CURPAR=',1P,11E10.3)
|
|---|
| 1005 | GOTO 7
|
|---|
| 1006 | ENDIF
|
|---|
| 1007 |
|
|---|
| 1008 | C GET NEXT PARTICLE FROM STACK, IF IRET=1 ALL PARTICLES ARE DONE
|
|---|
| 1009 | 4 CONTINUE
|
|---|
| 1010 | IRET1 = 0
|
|---|
| 1011 | CALL FSTACK
|
|---|
| 1012 | IF ( FMUADD ) THEN
|
|---|
| 1013 | IF ( (CURPAR(1) .EQ. 5 .OR. CURPAR(1) .EQ. 6)
|
|---|
| 1014 | * .AND. IRET1 .EQ. 0 .AND. .NOT. FMUORG ) THEN
|
|---|
| 1015 | DO J = 1,MAXLEN
|
|---|
| 1016 | AMUPAR(J) = CURPAR(J)
|
|---|
| 1017 | ENDDO
|
|---|
| 1018 | IF(DEBUG)WRITE(MDEBUG,*)'MAIN : MUON STORED IN AMUPAR'
|
|---|
| 1019 | FMUORG = .TRUE.
|
|---|
| 1020 | ENDIF
|
|---|
| 1021 | ENDIF
|
|---|
| 1022 |
|
|---|
| 1023 | IF ( IRET1 .EQ. 0 ) GOTO 7
|
|---|
| 1024 |
|
|---|
| 1025 | C-----------------------------------------------------------------------
|
|---|
| 1026 | C FINISH SHOWER AND PRINT INFORMATION
|
|---|
| 1027 | CALL OUTEND
|
|---|
| 1028 |
|
|---|
| 1029 |
|
|---|
| 1030 | * IF ( DEBUG ) WRITE(MDEBUG,442) NPARTO
|
|---|
| 1031 | *442 FORMAT(' MAIN : NPARTO='/(' ',10F10.0))
|
|---|
| 1032 |
|
|---|
| 1033 | IF ( FPRINT .OR. DEBUG ) THEN
|
|---|
| 1034 | IOBSLV = MIN( 5, NOBSLV )
|
|---|
| 1035 | WRITE(MONIOU,54) (K,K=1,IOBSLV)
|
|---|
| 1036 | 54 FORMAT (/' PARTICLES AT DETECTOR LEVEL :'/
|
|---|
| 1037 | * ' FOR LEVEL ', 5I13)
|
|---|
| 1038 | WRITE(MONIOU,55) (OBSLEV(K),K=1,IOBSLV)
|
|---|
| 1039 | 55 FORMAT ( ' HEIGHT IN CM ',1P, 5E13.3/)
|
|---|
| 1040 | WRITE(MONIOU,776) 'PROTONS ',(NPROTO(K),K=1,IOBSLV)
|
|---|
| 1041 | WRITE(MONIOU,776) 'ANTIPROTONS ',(NPROTB(K),K=1,IOBSLV)
|
|---|
| 1042 | WRITE(MONIOU,776) 'NEUTRONS ',(NNEUTR(K),K=1,IOBSLV)
|
|---|
| 1043 | WRITE(MONIOU,776) 'ANTINEUTRONS ',(NNEUTB(K),K=1,IOBSLV)
|
|---|
| 1044 | WRITE(MONIOU,776) 'PHOTONS ',(NPHOTO(K),K=1,IOBSLV)
|
|---|
| 1045 | WRITE(MONIOU,776) 'ELECTRONS ',(NELECT(K),K=1,IOBSLV)
|
|---|
| 1046 | WRITE(MONIOU,776) 'POSITRONS ',(NPOSIT(K),K=1,IOBSLV)
|
|---|
| 1047 | WRITE(MONIOU,776) 'NEUTRINOS ',(NNU (K),K=1,IOBSLV)
|
|---|
| 1048 | WRITE(MONIOU,776) 'MU - ',(NMUM (K),K=1,IOBSLV)
|
|---|
| 1049 | WRITE(MONIOU,776) 'MU + ',(NMUP (K),K=1,IOBSLV)
|
|---|
| 1050 | WRITE(MONIOU,776) 'PI 0 ',(NPI0 (K),K=1,IOBSLV)
|
|---|
| 1051 | WRITE(MONIOU,776) 'PI - ',(NPIM (K),K=1,IOBSLV)
|
|---|
| 1052 | WRITE(MONIOU,776) 'PI + ',(NPIP (K),K=1,IOBSLV)
|
|---|
| 1053 | WRITE(MONIOU,776) 'K0L ',(NK0L (K),K=1,IOBSLV)
|
|---|
| 1054 | WRITE(MONIOU,776) 'K0S ',(NK0S (K),K=1,IOBSLV)
|
|---|
| 1055 | WRITE(MONIOU,776) 'K - ',(NKMI (K),K=1,IOBSLV)
|
|---|
| 1056 | WRITE(MONIOU,776) 'K + ',(NKPL (K),K=1,IOBSLV)
|
|---|
| 1057 | WRITE(MONIOU,776) 'STR. BARYONS ',(NHYP (K),K=1,IOBSLV)
|
|---|
| 1058 | WRITE(MONIOU,776) 'DEUTERONS ',(NDEUT (K),K=1,IOBSLV)
|
|---|
| 1059 | WRITE(MONIOU,776) 'TRITONS ',(NTRIT (K),K=1,IOBSLV)
|
|---|
| 1060 | WRITE(MONIOU,776) 'ALPHAS ',(NALPHA(K),K=1,IOBSLV)
|
|---|
| 1061 | WRITE(MONIOU,776) 'OTHER PARTIC.',(NOTHER(K),K=1,IOBSLV)
|
|---|
| 1062 | WRITE(MONIOU,*)
|
|---|
| 1063 | WRITE(MONIOU,776) 'DECAYED MUONS',MUOND
|
|---|
| 1064 | 776 FORMAT(' NO OF ',A13, '= ',5F13.0)
|
|---|
| 1065 |
|
|---|
| 1066 | IF ( NOBSLV .GT. 5 ) THEN
|
|---|
| 1067 | IOBSLV = NOBSLV
|
|---|
| 1068 | WRITE(MONIOU,54) (K,K=6,IOBSLV)
|
|---|
| 1069 | WRITE(MONIOU,55) (OBSLEV(K),K=6,IOBSLV)
|
|---|
| 1070 | WRITE(MONIOU,776) 'PROTONS ',(NPROTO(K),K=6,IOBSLV)
|
|---|
| 1071 | WRITE(MONIOU,776) 'ANTIPROTONS ',(NPROTB(K),K=6,IOBSLV)
|
|---|
| 1072 | WRITE(MONIOU,776) 'NEUTRONS ',(NNEUTR(K),K=6,IOBSLV)
|
|---|
| 1073 | WRITE(MONIOU,776) 'ANTINEUTRONS ',(NNEUTB(K),K=6,IOBSLV)
|
|---|
| 1074 | WRITE(MONIOU,776) 'PHOTONS ',(NPHOTO(K),K=6,IOBSLV)
|
|---|
| 1075 | WRITE(MONIOU,776) 'ELECTRONS ',(NELECT(K),K=6,IOBSLV)
|
|---|
| 1076 | WRITE(MONIOU,776) 'POSITRONS ',(NPOSIT(K),K=6,IOBSLV)
|
|---|
| 1077 | WRITE(MONIOU,776) 'NEUTRINOS ',(NNU (K),K=6,IOBSLV)
|
|---|
| 1078 | WRITE(MONIOU,776) 'MU - ',(NMUM (K),K=6,IOBSLV)
|
|---|
| 1079 | WRITE(MONIOU,776) 'MU + ',(NMUP (K),K=6,IOBSLV)
|
|---|
| 1080 | WRITE(MONIOU,776) 'PI 0 ',(NPI0 (K),K=6,IOBSLV)
|
|---|
| 1081 | WRITE(MONIOU,776) 'PI - ',(NPIM (K),K=6,IOBSLV)
|
|---|
| 1082 | WRITE(MONIOU,776) 'PI + ',(NPIP (K),K=6,IOBSLV)
|
|---|
| 1083 | WRITE(MONIOU,776) 'K0L ',(NK0L (K),K=6,IOBSLV)
|
|---|
| 1084 | WRITE(MONIOU,776) 'K0S ',(NK0S (K),K=6,IOBSLV)
|
|---|
| 1085 | WRITE(MONIOU,776) 'K - ',(NKMI (K),K=6,IOBSLV)
|
|---|
| 1086 | WRITE(MONIOU,776) 'K + ',(NKPL (K),K=6,IOBSLV)
|
|---|
| 1087 | WRITE(MONIOU,776) 'STR. BARYONS ',(NHYP (K),K=6,IOBSLV)
|
|---|
| 1088 | WRITE(MONIOU,776) 'DEUTERONS ',(NDEUT (K),K=6,IOBSLV)
|
|---|
| 1089 | WRITE(MONIOU,776) 'TRITONS ',(NTRIT (K),K=6,IOBSLV)
|
|---|
| 1090 | WRITE(MONIOU,776) 'ALPHAS ',(NALPHA(K),K=6,IOBSLV)
|
|---|
| 1091 | WRITE(MONIOU,776) 'OTHER PARTIC.',(NOTHER(K),K=6,IOBSLV)
|
|---|
| 1092 | WRITE(MONIOU,*)
|
|---|
| 1093 | ENDIF
|
|---|
| 1094 | ENDIF
|
|---|
| 1095 |
|
|---|
| 1096 | C ADD UP FOR MEAN VALUES
|
|---|
| 1097 | DO 779 K = 1,25
|
|---|
| 1098 | DO 779 J = 1,10
|
|---|
| 1099 | MPARTO(J,K) = MPARTO(J,K) + NPARTO(J,K)
|
|---|
| 1100 | MPART2(J,K) = MPART2(J,K) + NPARTO(J,K)**2
|
|---|
| 1101 | 779 CONTINUE
|
|---|
| 1102 | EVTE(2) = SHOWNO
|
|---|
| 1103 | DO 335 K = 1,NOBSLV
|
|---|
| 1104 | EVTE(3) = EVTE(3) + NPHOTO(K)
|
|---|
| 1105 | EVTE(4) = EVTE(4) + NELECT(K) + NPOSIT(K)
|
|---|
| 1106 | EVTE(5) = EVTE(5) + NPROTO(K) + NPROTB(K) + NNEUTR(K) +
|
|---|
| 1107 | * NNEUTB(K) + NPI0(K) + NPIM(K) + NPIP(K) + NK0L(K) +
|
|---|
| 1108 | * NK0S(K) + NKMI(K) + NKPL(K) + NHYP(K) +
|
|---|
| 1109 | * NDEUT(K) + NTRIT(K) + NALPHA(K) + NOTHER(K)
|
|---|
| 1110 | EVTE(6) = EVTE(6) + NMUP(K) + NMUM(K)
|
|---|
| 1111 | 335 CONTINUE
|
|---|
| 1112 | EVTE(7) = NOPART
|
|---|
| 1113 |
|
|---|
| 1114 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,110)
|
|---|
| 1115 | * IFINNU,IFINPI,IFINET,IFINKA,IFINHY,
|
|---|
| 1116 | * IFINNU+IFINPI+IFINET+IFINKA+IFINHY,ELAST
|
|---|
| 1117 | 110 FORMAT(/' NO OF NUCLEONS PRODUCED IN FIRST INTERACTION =',I10/
|
|---|
| 1118 | * ' NO OF PIONS PRODUCED IN FIRST INTERACTION =',I10/
|
|---|
| 1119 | * ' NO OF ETAS PRODUCED IN FIRST INTERACTION =',I10/
|
|---|
| 1120 | * ' NO OF KAONS PRODUCED IN FIRST INTERACTION =',I10/
|
|---|
| 1121 | * ' NO OF S.BARYONS PRODUCED IN FIRST INTERACTION =',I10/
|
|---|
| 1122 | * ' TOTAL MULTIPLICITY OF FIRST INTERACTION =',I10/
|
|---|
| 1123 | * ' ELASTICITY OF FIRST INTERACTION =',F10.4)
|
|---|
| 1124 |
|
|---|
| 1125 | C PRINT OUT NKG RESULT FOR ONE SHOWER IF SELECTED
|
|---|
| 1126 | IF ( FNKG ) CALL AVAGE
|
|---|
| 1127 |
|
|---|
| 1128 | IF ( LLONGI ) THEN
|
|---|
| 1129 | C TREAT LONGITUDINAL DISTRIBUTIONS
|
|---|
| 1130 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1131 | c calculated here again, 'cos it's rewrite I dont know where
|
|---|
| 1132 | LPCT1 = INT( THICK0 * THSTPI )
|
|---|
| 1133 | LPCT2 = INT( (THICK0 + PRMPAR(3)*CHI) * THSTPI )
|
|---|
| 1134 | LPCT2 = MIN(NSTEP,LPCT2)
|
|---|
| 1135 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1136 | DO 980 J = LPCT1,NSTEP
|
|---|
| 1137 | C ADD ELECTRONS, POSITRONS, MUONS AND NUCLEI TO THE CHARGED PARTICLES
|
|---|
| 1138 | PLONG(J,7) = PLONG(J,7) + PLONG(J,2) + PLONG(J,3)
|
|---|
| 1139 | * + PLONG(J,4) + PLONG(J,5) + PLONG(J,8)
|
|---|
| 1140 | C ADD UP FOR MEAN VALUES OF LONGITUDINAL DISTRIBUTION
|
|---|
| 1141 | DO 979 K = 1,9
|
|---|
| 1142 | APLONG(J,K) = APLONG(J,K) + PLONG(J,K)
|
|---|
| 1143 | SPLONG(J,K) = SPLONG(J,K) + PLONG(J,K)**2
|
|---|
| 1144 | 979 CONTINUE
|
|---|
| 1145 | 980 CONTINUE
|
|---|
| 1146 |
|
|---|
| 1147 | C PRINT LONGITUDINAL DISTRIBUTIONS PER SHOWER
|
|---|
| 1148 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,910) THSTEP,
|
|---|
| 1149 | * 'GAMMAS','POSITRONS','ELECTRONS','MU-','MU+','HADRONS',
|
|---|
| 1150 | * 'CHARGED','NUCLEI','CERENKOV',
|
|---|
| 1151 | * (J*THSTEP,(PLONG(J,K),K=1,9),J=LPCT1,NSTEP)
|
|---|
| 1152 | 910 FORMAT(/' ---------- LONGITUDINAL DISTRIBUTION IN STEPS OF ',
|
|---|
| 1153 | * F5.0,' G/CM**2 ----------------'/
|
|---|
| 1154 | * ' DEPTH ',3A12,3A11,A12,A11,A12/
|
|---|
| 1155 | * (' ',F6.0,F13.0,2F12.0,3F11.0,F12.0,F11.0,1P,E12.5,0P) )
|
|---|
| 1156 | CJOK ADAPTED FOR HEAT CALCULATION
|
|---|
| 1157 | C910 FORMAT(/
|
|---|
| 1158 | C * ' LONGITUDINAL DISTRIBUTION IN STEPS OF ',F5.0,' G/CM**2'
|
|---|
| 1159 | C * /' ',92('=')/' DEPTH',8A10,A12/1P
|
|---|
| 1160 | C * (' ',0P,F6.0,1P,9E11.4))
|
|---|
| 1161 | CJOK
|
|---|
| 1162 |
|
|---|
| 1163 | IF ( FLGFIT ) THEN
|
|---|
| 1164 | C PERFORM FIT TO THE LONGITUDINAL DISTRIBUTION OF ALL CHARGED PARTICLES
|
|---|
| 1165 | C IF EGS IS SELECTED THIS IS THE DISTRIBUTION WHICH IS TO BE TAKEN
|
|---|
| 1166 | IF ( FEGS ) THEN
|
|---|
| 1167 | DO 930 J=0,NSTEP-LPCT1
|
|---|
| 1168 | DEP(J+1) = (J+LPCT1)*THSTEP
|
|---|
| 1169 | CHAPAR(J+1) = PLONG(J+LPCT1,7)
|
|---|
| 1170 | 930 CONTINUE
|
|---|
| 1171 | NSTP = NSTEP + 1 - LPCT1
|
|---|
| 1172 | WRITE(MONIOU,8229) 'ALL CHARGED PARTICLES'
|
|---|
| 1173 | 8229 FORMAT(/' FIT OF THE CURVE ',
|
|---|
| 1174 | * ' N(T) = P1*((T-P2)/(P3-P2))**((P3-T)/(P4+P5*T+P6*T**2))'/
|
|---|
| 1175 | * ' TO LONGITUDINAL DISTRIBUTION OF ',A35)
|
|---|
| 1176 | C IF NKG IS SELECTED ONLY THE ELECTRON DISTRIBUTION IS AVAILABLE
|
|---|
| 1177 | ELSEIF ( FNKG ) THEN
|
|---|
| 1178 | DEP(1) = 0.D0
|
|---|
| 1179 | CHAPAR(1) = 0.D0
|
|---|
| 1180 | DO 931 J = 1,IALT(1)
|
|---|
| 1181 | DEP(J+1) = TLEV(J)
|
|---|
| 1182 | CHAPAR(J+1) = SL(J)
|
|---|
| 1183 | 931 CONTINUE
|
|---|
| 1184 | NSTP = IALT(1) + 1
|
|---|
| 1185 | WRITE(MONIOU,8229) 'NKG ELECTRONS'
|
|---|
| 1186 | C IF NONE IS SELECTED IT DOES NOT REALLY MAKE SENSE TO FIT
|
|---|
| 1187 | C BUT LET'S TAKE THEN ALL CHARGED WHICH ARE MUONS AND HADRONS
|
|---|
| 1188 | ELSE
|
|---|
| 1189 | DO 932 J=0,NSTEP-LPCT1
|
|---|
| 1190 | DEP(J+1) = (J+LPCT1)*THSTEP
|
|---|
| 1191 | CHAPAR(J+1) = PLONG(J+LPCT1,7)
|
|---|
| 1192 | 932 CONTINUE
|
|---|
| 1193 | NSTP = NSTEP + 1 - LPCT1
|
|---|
| 1194 | WRITE(MONIOU,8229) 'MUONS AND CHARGED HADRONS'
|
|---|
| 1195 | ENDIF
|
|---|
| 1196 | IF ( NSTP .GT. 6 ) THEN
|
|---|
| 1197 | C THERE ARE MORE THAN 6 STEP VALUES, A FIT SHOULD BE POSSIBLE.
|
|---|
| 1198 | C DO THE FIT: NPAR AND FPARAM GIVE THE NUMBER OF PARAMETERS USED
|
|---|
| 1199 | C AND THE FINAL VALUES FOR THE PARAMETERS. CHISQ GIVES THE CHI**2/DOF
|
|---|
| 1200 | C FOR THE FIT.
|
|---|
| 1201 | CALL LONGFT(FPARAM,CHI2)
|
|---|
| 1202 | WRITE(MONIOU,8230) FPARAM,CHI2,CHI2/SQRT(FPARAM(1))*100.D0
|
|---|
| 1203 | 8230 FORMAT(' PARAMETERS = ',1P,6E12.4,0P/
|
|---|
| 1204 | * ' CHI**2/DOF = ',F10.1/
|
|---|
| 1205 | * ' AV. DEVIATION IN % = ',F10.4)
|
|---|
| 1206 | C STORE RESULT IN END EVENT BLOCK
|
|---|
| 1207 | DO 933 K = 1,6
|
|---|
| 1208 | EVTE(255+K) = FPARAM(K)
|
|---|
| 1209 | 933 CONTINUE
|
|---|
| 1210 | EVTE(262) = CHI2
|
|---|
| 1211 | ELSE
|
|---|
| 1212 | WRITE(MONIOU,*) 'NO LONGI. FIT POSSIBLE, ',
|
|---|
| 1213 | * ' NSTP = ',NSTP,' TOO SMALL.'
|
|---|
| 1214 | DO 934 K = 1,6
|
|---|
| 1215 | EVTE(255+K) = 0.
|
|---|
| 1216 | 934 CONTINUE
|
|---|
| 1217 | EVTE(262) = 0.
|
|---|
| 1218 | ENDIF
|
|---|
| 1219 | ENDIF
|
|---|
| 1220 | ENDIF
|
|---|
| 1221 |
|
|---|
| 1222 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1223 | c Saves statistics to sta###### file
|
|---|
| 1224 | call jcstadata(EVTH,EVTE,
|
|---|
| 1225 | + NPROTO,NPROTB,NNEUTR,NNEUTB,NPHOTO,NELECT,NPOSIT,
|
|---|
| 1226 | + NNU ,NMUM ,NMUP ,NPI0 ,NPIM ,NPIP ,NK0L ,
|
|---|
| 1227 | + NK0S ,NKMI ,NKPL ,NHYP ,NDEUT ,NTRIT ,NALPHA,
|
|---|
| 1228 | + NOTHER,IFINNU,IFINPI,IFINET,IFINKA,IFINHY,
|
|---|
| 1229 | + CERELE,CERHAD,PLONG,LPCT1,NSTEP,THSTEP)
|
|---|
| 1230 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1231 |
|
|---|
| 1232 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1233 | C WRITE SHOWER END TO OUTPUT BUFFER
|
|---|
| 1234 | c CALL TOBUF( EVTE,0 )
|
|---|
| 1235 | CALL TOBUF( EVTE,1 )
|
|---|
| 1236 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1237 | IF ( LCERFI ) THEN
|
|---|
| 1238 | CALL OUTND2
|
|---|
| 1239 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1240 | c CALL TOBUFC( EVTE,0 )
|
|---|
| 1241 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1242 | ENDIF
|
|---|
| 1243 |
|
|---|
| 1244 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,*)
|
|---|
| 1245 | * 'CERENKOV PH. FROM ELECTRONS = ',SNGL(CERELE),
|
|---|
| 1246 | * ' CERENKOV PH. FROM HADRONS = ',SNGL(CERHAD)
|
|---|
| 1247 | CERELE = 0.D0
|
|---|
| 1248 | CERHAD = 0.D0
|
|---|
| 1249 | NRECER = 0
|
|---|
| 1250 |
|
|---|
| 1251 | IF ( FPRINT .OR. DEBUG ) WRITE(MONIOU,210) SHOWNO
|
|---|
| 1252 | 210 FORMAT(/' END OF SHOWER NO ',I10)
|
|---|
| 1253 |
|
|---|
| 1254 | DO 19 J = 1,37
|
|---|
| 1255 | JNBIN(J) = JNBIN(J) + INBIN(J)
|
|---|
| 1256 | JPBIN(J) = JPBIN(J) + IPBIN(J)
|
|---|
| 1257 | JKBIN(J) = JKBIN(J) + IKBIN(J)
|
|---|
| 1258 | JHBIN(J) = JHBIN(J) + IHBIN(J)
|
|---|
| 1259 | 19 CONTINUE
|
|---|
| 1260 |
|
|---|
| 1261 | 2 CONTINUE
|
|---|
| 1262 | C END OF SHOWER LOOP
|
|---|
| 1263 |
|
|---|
| 1264 | C-----------------------------------------------------------------------
|
|---|
| 1265 | 992 CONTINUE
|
|---|
| 1266 |
|
|---|
| 1267 | C RESET NUMBER OF SHOWERS TO CORRECT VALUE
|
|---|
| 1268 | ISHW = I
|
|---|
| 1269 |
|
|---|
| 1270 | RUNE(3) = REAL(ISHW)
|
|---|
| 1271 | C WRITE RUN END TO OUTPUT BUFFER AND FINISH OUTPUT
|
|---|
| 1272 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1273 | c CALL TOBUF ( RUNE,1 )
|
|---|
| 1274 | call jcendrun(rune)
|
|---|
| 1275 | c IF ( LCERFI ) CALL TOBUFC( RUNE,1 )
|
|---|
| 1276 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1277 | C TIME SINCE BEGINNING
|
|---|
| 1278 | c ILEFTB = TIME()
|
|---|
| 1279 | ILEFTB = 1
|
|---|
| 1280 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1281 | TDIFF = ILEFTB - ILEFTA
|
|---|
| 1282 |
|
|---|
| 1283 | C MEAN VALUE FOR FIRST INTERACTION ALTITUDE (G/CM**2)
|
|---|
| 1284 | IF ( ISHW .GT. 1 ) THEN
|
|---|
| 1285 | CHISM2 = SQRT( ABS(CHISM2-CHISUM**2/ISHW) / (ISHW-1) )
|
|---|
| 1286 | CHISUM = CHISUM / ISHW
|
|---|
| 1287 | ELSE
|
|---|
| 1288 | CHISM2 = 0.D0
|
|---|
| 1289 | ENDIF
|
|---|
| 1290 |
|
|---|
| 1291 | C OUTPUTS FOR ALL SHOWERS
|
|---|
| 1292 | WRITE(MONIOU,201) ISHW,TDIFF,TDIFF/ISHW,IRECOR,IRECOR/ISHW,
|
|---|
| 1293 | * CHISUM,CHISM2
|
|---|
| 1294 | 201 FORMAT('1',10('='),' RUN SUMMARY ',56('=')//
|
|---|
| 1295 | * ' NUMBER OF GENERATED EVENTS = ',I10,/
|
|---|
| 1296 | * ' TOTAL TIME USED = ',E10.3,' SEC'/
|
|---|
| 1297 | * ' TIME PER EVENT = ',E10.3,' SEC'/
|
|---|
| 1298 | * ' TOTAL SPACE ON PATAPE USED = ',I10,' WORDS'/
|
|---|
| 1299 | * ' SPACE PER EVENT ON PATAPE = ',I10,' WORDS'/
|
|---|
| 1300 | * ' AVERAGE HEIGHT OF 1ST INT. = ',F10.3,' +-',F10.3,' G/CM**2'/)
|
|---|
| 1301 |
|
|---|
| 1302 | C ENERGY - MULTIPLICITY MATRIX FOR ALL SHOWERS
|
|---|
| 1303 | WRITE(MONIOU,209) (K,K=1,13),
|
|---|
| 1304 | * (J,(MULTOT(J,K),K=1,13),10**((J-4.)/3.),10**((J-3.)/3.),J=1,37),
|
|---|
| 1305 | * 1,(INT(10**((K-1.)/3.)+1),K = 2,13),
|
|---|
| 1306 | * 2,(INT(10**((K )/3.) ),K = 2,13)
|
|---|
| 1307 | 209 FORMAT(//' ENERGY - MULTIPLICITY MATRIX FOR ALL SHOWERS'/
|
|---|
| 1308 | * ' ENERGY RUNS VERTICALLY, MULTIPLICITY HORIZONTALLY'//,
|
|---|
| 1309 | * ' ',6X,5I9,3I8,5I7' ENERGY RANGE (GEV)'/
|
|---|
| 1310 | * 37(/' ',I4,1X,I10,4I9,3I8,5I7,1X,1P,2E10.1,0P)//
|
|---|
| 1311 | * ' MULT. ',5I9,3I8,5I7,4X,'LOWER BIN LIMIT'/
|
|---|
| 1312 | * ' RANGE ',5I9,3I8,5I7,4X,'UPPER BIN LIMIT')
|
|---|
| 1313 |
|
|---|
| 1314 | C GET MEAN OF ELASTICITY FOR ENERGY BINS
|
|---|
| 1315 | DO 3377 J = 1,37
|
|---|
| 1316 | NELMEA = 0
|
|---|
| 1317 | DO 3378 K = 1,10
|
|---|
| 1318 | NELMEA = NELMEA + IELDPA(J,K)
|
|---|
| 1319 | 3378 CONTINUE
|
|---|
| 1320 | IF ( NELMEA .NE. 0 ) ELMEAA(J) = ELMEAA(J) / NELMEA
|
|---|
| 1321 | 3377 CONTINUE
|
|---|
| 1322 |
|
|---|
| 1323 | C PRINT ENERGY - ELASTICITY MATRIX FOR ALL SHOWERS
|
|---|
| 1324 | WRITE(MONIOU,408) (K,K=1,10), (J,(IELDPA(J,K),K=1,10),
|
|---|
| 1325 | * ELMEAA(J),10**((J-4.D0)/3.D0),10**((J-3.)/3.D0),J=1,37),
|
|---|
| 1326 | * ((K-1)*0.1D0,K=1,10),(K*0.1D0,K=1,10)
|
|---|
| 1327 | 408 FORMAT (//' ENERGY - ELASTICITY MATRIX FOR ALL SHOWERS'/
|
|---|
| 1328 | * ' ENERGY RUNS VERTICALLY, ELASTICITY HORIZONTALLY'//
|
|---|
| 1329 | * ' ',5X,10I9,' MEAN EL. ENERGY RANGE (GEV)'/
|
|---|
| 1330 | * 37(/' ',I4,1X,10I9,2X,1P,E10.3,2E10.1,0P)//
|
|---|
| 1331 | * ' ELA. ',10F9.2,5X,'LOWER BIN LIMIT'/
|
|---|
| 1332 | * ' RANGE',10F9.2,5X,'UPPER BIN LIMIT')
|
|---|
| 1333 |
|
|---|
| 1334 | WRITE(MONIOU,204)
|
|---|
| 1335 | 204 FORMAT (//' INTERACTIONS PER KINETIC ENERGY INTERVAL FOR ALL ',
|
|---|
| 1336 | * 'SHOWERS'//' BIN LOWER LIMIT UPPER LIMIT ',
|
|---|
| 1337 | * 'NUCLEON PIONS KAONS S.BARYONS TOTAL'/
|
|---|
| 1338 | * 12X,'IN GEV',9X,'IN GEV',7X,
|
|---|
| 1339 | * ' EVENTS EVENTS EVENTS EVENTS '//)
|
|---|
| 1340 | WRITE(MONIOU,207) (I,SABIN(I),SBBIN(I),JNBIN(I),JPBIN(I),JKBIN(I)
|
|---|
| 1341 | * ,JHBIN(I),JNBIN(I)+JPBIN(I)+JKBIN(I)+JHBIN(I),I=1,37)
|
|---|
| 1342 | 207 FORMAT(' ',I5,1P,2E15.4,0P,I12,3I10,I11)
|
|---|
| 1343 |
|
|---|
| 1344 | IF ( .NOT.GHEISH ) THEN
|
|---|
| 1345 | C PRINT ELASTICITY STATISTICS
|
|---|
| 1346 | WRITE(MONIOU,89) (I,(I-1)*.05,I*.05,
|
|---|
| 1347 | * IELIS(I),IELHM(I),IELNU(I),IELPI(I),I = 1,20)
|
|---|
| 1348 | 89 FORMAT (//' ELASTICITY STATISTICS '//
|
|---|
| 1349 | * ' BIN LOW HIGH EDGE FOR ISOBARS HEAVY MESONS',
|
|---|
| 1350 | * ' SINGLE NUCLEONS AND PIONS'/
|
|---|
| 1351 | * (' ',I3,' ',F4.2,' ',F4.2,' ',4I17))
|
|---|
| 1352 | ENDIF
|
|---|
| 1353 |
|
|---|
| 1354 | C CALCULATE MEAN VALUES AND STANDARD DEVIATIONS OF PARTICLE NUMBERS
|
|---|
| 1355 | IF ( ISHW .GT. 1 ) THEN
|
|---|
| 1356 | DO 879 K = 1,25
|
|---|
| 1357 | DO 879 J = 1,NOBSLV
|
|---|
| 1358 | MPART2(J,K) = SQRT( abs(MPART2(J,K)-MPARTO(J,K)**2/ISHW)
|
|---|
| 1359 | * /(ISHW-1) )
|
|---|
| 1360 | MPARTO(J,K) = MPARTO(J,K)/ISHW
|
|---|
| 1361 | 879 CONTINUE
|
|---|
| 1362 | ELSE
|
|---|
| 1363 | DO 880 K = 1,25
|
|---|
| 1364 | DO 880 J = 1,NOBSLV
|
|---|
| 1365 | MPART2(J,K) = 0.D0
|
|---|
| 1366 | 880 CONTINUE
|
|---|
| 1367 | ENDIF
|
|---|
| 1368 |
|
|---|
| 1369 | C PRINT MEAN VALUES AND STANDARD DEVIATIONS OF PARTICLE NUMBERS
|
|---|
| 1370 | IOBSLV = MIN( 3, NOBSLV )
|
|---|
| 1371 | WRITE(MONIOU,854) (K,K=1,IOBSLV)
|
|---|
| 1372 | 854 FORMAT (/ ' AVERAGE NUMBER OF PARTICLES PER EVENT :'/
|
|---|
| 1373 | * ' FROM LEVEL NUMBER ', 3(10X,I10,10X) )
|
|---|
| 1374 | WRITE(MONIOU,855) (OBSLEV(K),K=1,IOBSLV)
|
|---|
| 1375 | 855 FORMAT ( ' HEIGHT IN CM',1P,3(20X,E10.3)/)
|
|---|
| 1376 |
|
|---|
| 1377 | WRITE(MONIOU,778)'PROTONS ',(MPROTO(K),MPROT2(K),K=1,IOBSLV)
|
|---|
| 1378 | WRITE(MONIOU,778)'ANTIPROTONS ',(MPROTB(K),MPRTB2(K),K=1,IOBSLV)
|
|---|
| 1379 | WRITE(MONIOU,778)'NEUTRONS ',(MNEUTR(K),MNETR2(K),K=1,IOBSLV)
|
|---|
| 1380 | WRITE(MONIOU,778)'ANTINEUTRONS',(MNEUTB(K),MNETB2(K),K=1,IOBSLV)
|
|---|
| 1381 | WRITE(MONIOU,778)'PHOTONS ',(MPHOTO(K),MPHOT2(K),K=1,IOBSLV)
|
|---|
| 1382 | WRITE(MONIOU,778)'ELECTRONS ',(MELECT(K),MELEC2(K),K=1,IOBSLV)
|
|---|
| 1383 | WRITE(MONIOU,778)'POSITRONS ',(MPOSIT(K),MPOSI2(K),K=1,IOBSLV)
|
|---|
| 1384 | WRITE(MONIOU,778)'NEUTRINOS ',(MNU (K),MNU2 (K),K=1,IOBSLV)
|
|---|
| 1385 | WRITE(MONIOU,778)'MU - ',(MMUM (K),MMUM2 (K),K=1,IOBSLV)
|
|---|
| 1386 | WRITE(MONIOU,778)'MU + ',(MMUP (K),MMUP2 (K),K=1,IOBSLV)
|
|---|
| 1387 | WRITE(MONIOU,778)'PI 0 ',(MPI0 (K),MPI02 (K),K=1,IOBSLV)
|
|---|
| 1388 | WRITE(MONIOU,778)'PI - ',(MPIM (K),MPIM2 (K),K=1,IOBSLV)
|
|---|
| 1389 | WRITE(MONIOU,778)'PI + ',(MPIP (K),MPIP2 (K),K=1,IOBSLV)
|
|---|
| 1390 | WRITE(MONIOU,778)'K0L ',(MK0L (K),MK0L2 (K),K=1,IOBSLV)
|
|---|
| 1391 | WRITE(MONIOU,778)'K0S ',(MK0S (K),MK0S2 (K),K=1,IOBSLV)
|
|---|
| 1392 | WRITE(MONIOU,778)'K - ',(MKMI (K),MKMI2 (K),K=1,IOBSLV)
|
|---|
| 1393 | WRITE(MONIOU,778)'K + ',(MKPL (K),MKPL2 (K),K=1,IOBSLV)
|
|---|
| 1394 | WRITE(MONIOU,778)'STR. BARYONS',(MHYP (K),MHYP2 (K),K=1,IOBSLV)
|
|---|
| 1395 | WRITE(MONIOU,778)'DEUTERONS ',(MDEUT (K),MDEUT2(K),K=1,IOBSLV)
|
|---|
| 1396 | WRITE(MONIOU,778)'TRITONS ',(MTRIT (K),MTRIT2(K),K=1,IOBSLV)
|
|---|
| 1397 | WRITE(MONIOU,778)'ALPHAS ',(MALPHA(K),MALPH2(K),K=1,IOBSLV)
|
|---|
| 1398 | WRITE(MONIOU,778)'OTHER PART. ',(MOTHER(K),MOTH2 (K),K=1,IOBSLV)
|
|---|
| 1399 | WRITE(MONIOU,*)
|
|---|
| 1400 | 778 FORMAT(' NO OF ',A12,' = ',3(F13.1,' +-',F13.1,' '))
|
|---|
| 1401 |
|
|---|
| 1402 | IF ( NOBSLV .GT. 3 ) THEN
|
|---|
| 1403 | IOBSLV = MIN( 6, NOBSLV )
|
|---|
| 1404 | WRITE(MONIOU,854) (K,K=4,IOBSLV)
|
|---|
| 1405 | WRITE(MONIOU,855) (OBSLEV(K),K=4,IOBSLV)
|
|---|
| 1406 | WRITE(MONIOU,778)'PROTONS ',(MPROTO(K),MPROT2(K),K=4,IOBSLV)
|
|---|
| 1407 | WRITE(MONIOU,778)'ANTIPROTONS ',(MPROTB(K),MPRTB2(K),K=4,IOBSLV)
|
|---|
| 1408 | WRITE(MONIOU,778)'NEUTRONS ',(MNEUTR(K),MNETR2(K),K=4,IOBSLV)
|
|---|
| 1409 | WRITE(MONIOU,778)'ANTINEUTRONS',(MNEUTB(K),MNETB2(K),K=4,IOBSLV)
|
|---|
| 1410 | WRITE(MONIOU,778)'PHOTONS ',(MPHOTO(K),MPHOT2(K),K=4,IOBSLV)
|
|---|
| 1411 | WRITE(MONIOU,778)'ELECTRONS ',(MELECT(K),MELEC2(K),K=4,IOBSLV)
|
|---|
| 1412 | WRITE(MONIOU,778)'POSITRONS ',(MPOSIT(K),MPOSI2(K),K=4,IOBSLV)
|
|---|
| 1413 | WRITE(MONIOU,778)'NEUTRINOS ',(MNU (K),MNU2 (K),K=4,IOBSLV)
|
|---|
| 1414 | WRITE(MONIOU,778)'MU - ',(MMUM (K),MMUM2 (K),K=4,IOBSLV)
|
|---|
| 1415 | WRITE(MONIOU,778)'MU + ',(MMUP (K),MMUP2 (K),K=4,IOBSLV)
|
|---|
| 1416 | WRITE(MONIOU,778)'PI 0 ',(MPI0 (K),MPI02 (K),K=4,IOBSLV)
|
|---|
| 1417 | WRITE(MONIOU,778)'PI - ',(MPIM (K),MPIM2 (K),K=4,IOBSLV)
|
|---|
| 1418 | WRITE(MONIOU,778)'PI + ',(MPIP (K),MPIP2 (K),K=4,IOBSLV)
|
|---|
| 1419 | WRITE(MONIOU,778)'K0L ',(MK0L (K),MK0L2 (K),K=4,IOBSLV)
|
|---|
| 1420 | WRITE(MONIOU,778)'K0S ',(MK0S (K),MK0S2 (K),K=4,IOBSLV)
|
|---|
| 1421 | WRITE(MONIOU,778)'K - ',(MKMI (K),MKMI2 (K),K=4,IOBSLV)
|
|---|
| 1422 | WRITE(MONIOU,778)'K + ',(MKPL (K),MKPL2 (K),K=4,IOBSLV)
|
|---|
| 1423 | WRITE(MONIOU,778)'STR. BARYONS',(MHYP (K),MHYP2 (K),K=4,IOBSLV)
|
|---|
| 1424 | WRITE(MONIOU,778)'DEUTERONS ',(MDEUT (K),MDEUT2(K),K=4,IOBSLV)
|
|---|
| 1425 | WRITE(MONIOU,778)'TRITONS ',(MTRIT (K),MTRIT2(K),K=4,IOBSLV)
|
|---|
| 1426 | WRITE(MONIOU,778)'ALPHAS ',(MALPHA(K),MALPH2(K),K=4,IOBSLV)
|
|---|
| 1427 | WRITE(MONIOU,778)'OTHER PART. ',(MOTHER(K),MOTH2 (K),K=4,IOBSLV)
|
|---|
| 1428 | WRITE(MONIOU,*)
|
|---|
| 1429 |
|
|---|
| 1430 | IF ( NOBSLV .GT. 6 ) THEN
|
|---|
| 1431 | IOBSLV = MIN( 9, NOBSLV )
|
|---|
| 1432 | WRITE(MONIOU,854) (K,K=7,IOBSLV)
|
|---|
| 1433 | WRITE(MONIOU,855) (OBSLEV(K),K=7,IOBSLV)
|
|---|
| 1434 | WRITE(MONIOU,778)'PROTONS ',(MPROTO(K),MPROT2(K),K=7,IOBSLV)
|
|---|
| 1435 | WRITE(MONIOU,778)'ANTIPROTONS ',(MPROTB(K),MPRTB2(K),K=7,IOBSLV)
|
|---|
| 1436 | WRITE(MONIOU,778)'NEUTRONS ',(MNEUTR(K),MNETR2(K),K=7,IOBSLV)
|
|---|
| 1437 | WRITE(MONIOU,778)'ANTINEUTRONS',(MNEUTB(K),MNETB2(K),K=7,IOBSLV)
|
|---|
| 1438 | WRITE(MONIOU,778)'PHOTONS ',(MPHOTO(K),MPHOT2(K),K=7,IOBSLV)
|
|---|
| 1439 | WRITE(MONIOU,778)'ELECTRONS ',(MELECT(K),MELEC2(K),K=7,IOBSLV)
|
|---|
| 1440 | WRITE(MONIOU,778)'POSITRONS ',(MPOSIT(K),MPOSI2(K),K=7,IOBSLV)
|
|---|
| 1441 | WRITE(MONIOU,778)'NEUTRINOS ',(MNU (K),MNU2 (K),K=7,IOBSLV)
|
|---|
| 1442 | WRITE(MONIOU,778)'MU - ',(MMUM (K),MMUM2 (K),K=7,IOBSLV)
|
|---|
| 1443 | WRITE(MONIOU,778)'MU + ',(MMUP (K),MMUP2 (K),K=7,IOBSLV)
|
|---|
| 1444 | WRITE(MONIOU,778)'PI 0 ',(MPI0 (K),MPI02 (K),K=7,IOBSLV)
|
|---|
| 1445 | WRITE(MONIOU,778)'PI - ',(MPIM (K),MPIM2 (K),K=7,IOBSLV)
|
|---|
| 1446 | WRITE(MONIOU,778)'PI + ',(MPIP (K),MPIP2 (K),K=7,IOBSLV)
|
|---|
| 1447 | WRITE(MONIOU,778)'K0L ',(MK0L (K),MK0L2 (K),K=7,IOBSLV)
|
|---|
| 1448 | WRITE(MONIOU,778)'K0S ',(MK0S (K),MK0S2 (K),K=7,IOBSLV)
|
|---|
| 1449 | WRITE(MONIOU,778)'K - ',(MKMI (K),MKMI2 (K),K=7,IOBSLV)
|
|---|
| 1450 | WRITE(MONIOU,778)'K + ',(MKPL (K),MKPL2 (K),K=7,IOBSLV)
|
|---|
| 1451 | WRITE(MONIOU,778)'STR. BARYONS',(MHYP (K),MHYP2 (K),K=7,IOBSLV)
|
|---|
| 1452 | WRITE(MONIOU,778)'DEUTERONS ',(MDEUT (K),MDEUT2(K),K=7,IOBSLV)
|
|---|
| 1453 | WRITE(MONIOU,778)'TRITONS ',(MTRIT (K),MTRIT2(K),K=7,IOBSLV)
|
|---|
| 1454 | WRITE(MONIOU,778)'ALPHAS ',(MALPHA(K),MALPH2(K),K=7,IOBSLV)
|
|---|
| 1455 | WRITE(MONIOU,778)'OTHER PART. ',(MOTHER(K),MOTH2 (K),K=7,IOBSLV)
|
|---|
| 1456 | WRITE(MONIOU,*)
|
|---|
| 1457 |
|
|---|
| 1458 | IF ( NOBSLV .GT. 9 ) THEN
|
|---|
| 1459 | IOBSLV = MIN( 10, NOBSLV )
|
|---|
| 1460 | WRITE(MONIOU,854) (K,K=9,IOBSLV)
|
|---|
| 1461 | WRITE(MONIOU,855) (OBSLEV(K),K=9,IOBSLV)
|
|---|
| 1462 | WRITE(MONIOU,778)'PROTONS ',(MPROTO(K),MPROT2(K),K=9,IOBSLV)
|
|---|
| 1463 | WRITE(MONIOU,778)'ANTIPROTONS ',(MPROTB(K),MPRTB2(K),K=9,IOBSLV)
|
|---|
| 1464 | WRITE(MONIOU,778)'NEUTRONS ',(MNEUTR(K),MNETR2(K),K=9,IOBSLV)
|
|---|
| 1465 | WRITE(MONIOU,778)'ANTINEUTRONS',(MNEUTB(K),MNETB2(K),K=9,IOBSLV)
|
|---|
| 1466 | WRITE(MONIOU,778)'PHOTONS ',(MPHOTO(K),MPHOT2(K),K=9,IOBSLV)
|
|---|
| 1467 | WRITE(MONIOU,778)'ELECTRONS ',(MELECT(K),MELEC2(K),K=9,IOBSLV)
|
|---|
| 1468 | WRITE(MONIOU,778)'POSITRONS ',(MPOSIT(K),MPOSI2(K),K=9,IOBSLV)
|
|---|
| 1469 | WRITE(MONIOU,778)'NEUTRINOS ',(MNU (K),MNU2 (K),K=9,IOBSLV)
|
|---|
| 1470 | WRITE(MONIOU,778)'MU - ',(MMUM (K),MMUM2 (K),K=9,IOBSLV)
|
|---|
| 1471 | WRITE(MONIOU,778)'MU + ',(MMUP (K),MMUP2 (K),K=9,IOBSLV)
|
|---|
| 1472 | WRITE(MONIOU,778)'PI 0 ',(MPI0 (K),MPI02 (K),K=9,IOBSLV)
|
|---|
| 1473 | WRITE(MONIOU,778)'PI - ',(MPIM (K),MPIM2 (K),K=9,IOBSLV)
|
|---|
| 1474 | WRITE(MONIOU,778)'PI + ',(MPIP (K),MPIP2 (K),K=9,IOBSLV)
|
|---|
| 1475 | WRITE(MONIOU,778)'K0L ',(MK0L (K),MK0L2 (K),K=9,IOBSLV)
|
|---|
| 1476 | WRITE(MONIOU,778)'K0S ',(MK0S (K),MK0S2 (K),K=9,IOBSLV)
|
|---|
| 1477 | WRITE(MONIOU,778)'K - ',(MKMI (K),MKMI2 (K),K=9,IOBSLV)
|
|---|
| 1478 | WRITE(MONIOU,778)'K + ',(MKPL (K),MKPL2 (K),K=9,IOBSLV)
|
|---|
| 1479 | WRITE(MONIOU,778)'STR. BARYONS',(MHYP (K),MHYP2 (K),K=9,IOBSLV)
|
|---|
| 1480 | WRITE(MONIOU,778)'DEUTERONS ',(MDEUT (K),MDEUT2(K),K=9,IOBSLV)
|
|---|
| 1481 | WRITE(MONIOU,778)'TRITONS ',(MTRIT (K),MTRIT2(K),K=9,IOBSLV)
|
|---|
| 1482 | WRITE(MONIOU,778)'ALPHAS ',(MALPHA(K),MALPH2(K),K=9,IOBSLV)
|
|---|
| 1483 | WRITE(MONIOU,778)'OTHER PART. ',(MOTHER(K),MOTH2 (K),K=9,IOBSLV)
|
|---|
| 1484 | WRITE(MONIOU,*)
|
|---|
| 1485 | ENDIF
|
|---|
| 1486 |
|
|---|
| 1487 | ENDIF
|
|---|
| 1488 | ENDIF
|
|---|
| 1489 |
|
|---|
| 1490 | C PRINT OUT NKG RESULT FOR ALL SHOWERS IF SELECTED
|
|---|
| 1491 | IF ( FNKG ) CALL MITAGE
|
|---|
| 1492 |
|
|---|
| 1493 | C CALCULATE MEAN VALUES AND SIGMAS OF LONGITUDINAL DISTRIBUTION
|
|---|
| 1494 | IF ( LLONGI ) THEN
|
|---|
| 1495 | IF ( ISHW .GT. 1 ) THEN
|
|---|
| 1496 | DO 790 K = 1,9
|
|---|
| 1497 | DO 789 J = LPCT1,NSTEP
|
|---|
| 1498 | SPLONG(J,K) = SQRT( abs(SPLONG(J,K)-APLONG(J,K)**2/ISHW)
|
|---|
| 1499 | * /(ISHW-1) )
|
|---|
| 1500 | APLONG(J,K) = APLONG(J,K)/ISHW
|
|---|
| 1501 | 789 CONTINUE
|
|---|
| 1502 | 790 CONTINUE
|
|---|
| 1503 | ELSE
|
|---|
| 1504 | DO 990 K = 1,9
|
|---|
| 1505 | DO 989 J = LPCT1,NSTEP
|
|---|
| 1506 | SPLONG(J,K) = 0.D0
|
|---|
| 1507 | 989 CONTINUE
|
|---|
| 1508 | 990 CONTINUE
|
|---|
| 1509 | ENDIF
|
|---|
| 1510 |
|
|---|
| 1511 | C PRINT AVERAGE LONGITUDINAL DISTRIBUTIONS
|
|---|
| 1512 | WRITE(MONIOU,911) THSTEP,
|
|---|
| 1513 | * 'GAMMAS ','POSITRONS','ELECTRONS','MU- ','MU+ ',
|
|---|
| 1514 | * (J*THSTEP,(APLONG(J,K),SPLONG(J,K),K=1,5),J=LPCT1,NSTEP)
|
|---|
| 1515 | 911 FORMAT(/' AVERAGE LONGITUDINAL DISTRIBUTION IN STEPS OF ',
|
|---|
| 1516 | * F5.0,' G/CM**2 '/' ',131('=')/
|
|---|
| 1517 | * ' DEPTH',8X,3(A10,16X),A9,15X,A9 //
|
|---|
| 1518 | * (' ',F5.0,2X,1P,E11.4,'+-',E11.4,0P,1X,F12.0,'+-',F11.0,
|
|---|
| 1519 | * 1X,F12.0,'+-',E11.4,1X,F11.1,'+-',F10.1,
|
|---|
| 1520 | * 1X,F11.1,'+-',F10.1 ))
|
|---|
| 1521 | WRITE(MONIOU,912) THSTEP,
|
|---|
| 1522 | * 'HADRONS','CHARGED','NUCLEI','CERENKOV',
|
|---|
| 1523 | * (J*THSTEP,(APLONG(J,K),SPLONG(J,K),K=6,9),J=LPCT1,NSTEP)
|
|---|
| 1524 | 912 FORMAT(/' AVERAGE LONGITUDINAL DISTRIBUTION IN STEPS OF ',
|
|---|
| 1525 | * F5.0,' G/CM**2 '/' ',115('=')/
|
|---|
| 1526 | * ' DEPTH',8X,A9,16X,A10,16X,A9,21X,A9 //
|
|---|
| 1527 | * (' ',F5.0,1X,F11.1,'+-',F11.1,1X,F12.0,'+-',F12.0,
|
|---|
| 1528 | * 2X,F10.1,'+-',F10.1,1X,1P,E16.6,'+-',E16.6,0P))
|
|---|
| 1529 | ENDIF
|
|---|
| 1530 |
|
|---|
| 1531 | IF ( FLGFIT ) THEN
|
|---|
| 1532 | C PERFORM FIT TO THE LONGITUDINAL DISTRIBUTION OF ALL CHARGED PARTICLES
|
|---|
| 1533 | C IF EGS IS SELECTED THIS IS THE DISTRIBUTION WHICH IS TO BE TAKEN
|
|---|
| 1534 | IF ( FEGS ) THEN
|
|---|
| 1535 | DO 730 J=0,NSTEP-LPCT1
|
|---|
| 1536 | DEP(J+1) = (J+LPCT1)*THSTEP
|
|---|
| 1537 | CHAPAR(J+1) = APLONG(J+LPCT1,7)
|
|---|
| 1538 | 730 CONTINUE
|
|---|
| 1539 | NSTP = NSTEP + 1 - LPCT1
|
|---|
| 1540 | WRITE(MONIOU,8229) 'AVERAGE ALL CHARGED PARTICLES'
|
|---|
| 1541 | C IF NKG IS SELECTED ONLY THE ELECTRON DISTRIBUTION IS AVAILABLE
|
|---|
| 1542 | ELSEIF ( FNKG ) THEN
|
|---|
| 1543 | DEP(1) = 0.D0
|
|---|
| 1544 | CHAPAR(1) = 0.D0
|
|---|
| 1545 | DO 731 J = 1,IALT(1)
|
|---|
| 1546 | DEP(J+1) = TLEV(J)
|
|---|
| 1547 | CHAPAR(J+1) = SEL(J)/ISHW
|
|---|
| 1548 | 731 CONTINUE
|
|---|
| 1549 | NSTP = IALT(1) + 1
|
|---|
| 1550 | WRITE(MONIOU,8229) 'AVERAGE NKG ELECTRONS'
|
|---|
| 1551 | C IF NONE IS SELECTED IT DOES NOT REALLY MAKE SENSE TO FIT
|
|---|
| 1552 | C BUT LET'S TAKE THEN ALL CHARGED WHICH ARE MUONS AND HADRONS
|
|---|
| 1553 | ELSE
|
|---|
| 1554 | DO 732 J=0,NSTEP-LPCT1
|
|---|
| 1555 | DEP(J+1) = (J+LPCT1)*THSTEP
|
|---|
| 1556 | CHAPAR(J+1) = APLONG(J+LPCT1,7)
|
|---|
| 1557 | 732 CONTINUE
|
|---|
| 1558 | NSTP = NSTEP + 1 - LPCT1
|
|---|
| 1559 | WRITE(MONIOU,8229) 'AVERAGE MUONS AND CHARGED HADRONS'
|
|---|
| 1560 | ENDIF
|
|---|
| 1561 | IF ( NSTP .GT. 6 ) THEN
|
|---|
| 1562 | C THERE ARE MORE THAN 6 STEP VALUES, A FIT SHOULD BE POSSIBLE.
|
|---|
| 1563 | C DO THE FIT: NPAR AND FPARAM GIVE THE NUMBER OF PARAMETERS USED
|
|---|
| 1564 | C AND THE FINAL VALUES FOR THE PARAMETERS. CHISQ GIVES THE CHI**2/DOF
|
|---|
| 1565 | C FOR THE FIT.
|
|---|
| 1566 | CALL LONGFT(FPARAM,CHI2)
|
|---|
| 1567 | WRITE(MONIOU,8230) FPARAM,CHI2,CHI2/SQRT(FPARAM(1))*100.D0
|
|---|
| 1568 | ELSE
|
|---|
| 1569 | WRITE(MONIOU,*) 'NO LONGI. FIT POSSIBLE, ',
|
|---|
| 1570 | * ' NSTP = ',NSTP,' TOO SMALL.'
|
|---|
| 1571 | ENDIF
|
|---|
| 1572 | ENDIF
|
|---|
| 1573 |
|
|---|
| 1574 |
|
|---|
| 1575 | C CONTROL PRINT OUTPUT OF CONSTANTS
|
|---|
| 1576 | IF ( DEBUG ) THEN
|
|---|
| 1577 | CALL STAEND
|
|---|
| 1578 | WRITE(MDEBUG,*) 'MAIN : STAEND CALLED'
|
|---|
| 1579 | ENDIF
|
|---|
| 1580 |
|
|---|
| 1581 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1582 | call jcenddata(runh,rune)
|
|---|
| 1583 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1584 |
|
|---|
| 1585 | WRITE(MONIOU,*)' '
|
|---|
| 1586 | CALL PRTIME(TTIME)
|
|---|
| 1587 | WRITE(MONIOU,101)
|
|---|
| 1588 | 101 FORMAT (/' ',10('='),' END OF RUN ',67('='))
|
|---|
| 1589 |
|
|---|
| 1590 | C CLOSE ALL OPEN UNITS
|
|---|
| 1591 | IF ( MONIOU .NE. 6 ) CLOSE( MONIOU )
|
|---|
| 1592 | IF ( MDEBUG .NE. 6 ) CLOSE( MDEBUG )
|
|---|
| 1593 | CLOSE( EXST )
|
|---|
| 1594 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1595 | c CLOSE( PATAPE )
|
|---|
| 1596 | c IF ( LCERFI ) CLOSE( CETAPE )
|
|---|
| 1597 | c>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
|
|---|
| 1598 |
|
|---|
| 1599 | STOP
|
|---|
| 1600 | END
|
|---|
