source: trunk/MagicSoft/Simulation/Corsika/Mmcs/start.f@ 6724

Last change on this file since 6724 was 287, checked in by harald, 25 years ago
Some minor changes in the files main.f and start.f due to some problems during the compilation on a linux computer (pchegra4.mppmu.mpg.de). After this changes it is possible to compile the program under linux, but there are still a lot of warnings! A small discussion with Jose Carlos and Aitor leads to the conclusion that the executable is good for mass production!!
File size: 24.5 KB
Line 
1 SUBROUTINE START
2
3C-----------------------------------------------------------------------
4C START
5C
6C PERFORMS INITIALISATIONS AND CHECKS AT THE BEGINNING OF RUN.
7C CALLS DATAC TO READ IN DATA CARDS.
8C CHECKS AND INITIALIZES SELECTED HADRONIC INTERACTION MODEL.
9C THIS SUBROUTINE IS CALLED FROM MAIN
10C
11C REDESIGN: J. KNAPP IK1 FZK KARLSRUHE
12C-----------------------------------------------------------------------
13
14 IMPLICIT NONE
15*KEEP,AIR.
16 COMMON /AIR/ COMPOS,PROBTA,AVERAW,AVOGAD
17 DOUBLE PRECISION COMPOS(3),PROBTA(3),AVERAW,AVOGAD
18*KEEP,ANNI.
19 COMMON /ANNI/ CAN,CANN
20 DOUBLE PRECISION CAN(50),CANN(50)
21*KEEP,ATMOS.
22 COMMON /ATMOS/ AATM,BATM,CATM,DATM
23 DOUBLE PRECISION AATM(5),BATM(5),CATM(5),DATM(5)
24*KEEP,ATMOS2.
25 COMMON /ATMOS2/ HLAY,THICKL
26 DOUBLE PRECISION HLAY(5),THICKL(5)
27*KEEP,BUFFS.
28 COMMON /BUFFS/ RUNH,RUNE,EVTH,EVTE,DATAB,LH
29 INTEGER MAXBUF,MAXLEN
30 PARAMETER (MAXBUF=39*7)
31 PARAMETER (MAXLEN=12)
32 REAL RUNH(MAXBUF),EVTH(MAXBUF),EVTE(MAXBUF),
33 * RUNE(MAXBUF),DATAB(MAXBUF)
34 INTEGER LH
35 CHARACTER*4 CRUNH,CRUNE,CEVTH,CEVTE
36 EQUIVALENCE (RUNH(1),CRUNH), (RUNE(1),CRUNE)
37 EQUIVALENCE (EVTH(1),CEVTH), (EVTE(1),CEVTE)
38*KEEP,CONST.
39 COMMON /CONST/ PI,PI2,OB3,TB3,ENEPER
40 DOUBLE PRECISION PI,PI2,OB3,TB3,ENEPER
41*KEEP,DPMFLG.
42 COMMON /DPMFLG/ NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
43 INTEGER NFLAIN,NFLDIF,NFLPI0,NFLCHE,NFLPIF,NFRAGM
44*KEEP,EDECAY.
45 COMMON /EDECAY/ CETA
46 DOUBLE PRECISION CETA(5)
47*KEEP,ELABCT.
48 COMMON /ELABCT/ ELCUT
49 DOUBLE PRECISION ELCUT(4)
50*KEEP,ETHMAP.
51 COMMON /ETHMAP/ ECTMAP,ELEFT
52 DOUBLE PRECISION ECTMAP,ELEFT
53*KEEP,KAONS.
54 COMMON /KAONS/ CKA
55 DOUBLE PRECISION CKA(80)
56*KEEP,MAGNET.
57 COMMON /MAGNET/ BX,BZ,BVAL,BNORMC,BNORM,COSB,SINB,BLIMIT
58 DOUBLE PRECISION BX,BZ,BVAL,BNORMC
59 REAL BNORM,COSB,SINB,BLIMIT
60*KEEP,MUMULT.
61 COMMON /MUMULT/ CHC,OMC,FMOLI
62 DOUBLE PRECISION CHC,OMC
63 LOGICAL FMOLI
64*KEEP,MUPART.
65 COMMON /MUPART/ AMUPAR,BCUT,CMUON,FMUBRM,FMUORG
66 DOUBLE PRECISION AMUPAR(14),BCUT,CMUON(11)
67 LOGICAL FMUBRM,FMUORG
68*KEEP,NCSNCS.
69 COMMON /NCSNCS/ SIGN30,SIGN45,SIGN60,SIGO30,SIGO45,SIGO60,
70 * SIGA30,SIGA45,SIGA60,PNOA30,PNOA45,PNOA60,
71 * SIG30A,SIG45A,SIG60A
72 DOUBLE PRECISION SIGN30(56),SIGN45(56),SIGN60(56),
73 * SIGO30(56),SIGO45(56),SIGO60(56),
74 * SIGA30(56),SIGA45(56),SIGA60(56),
75 * PNOA30(1540,3),PNOA45(1540,3),PNOA60(1540,3),
76 * SIG30A(56),SIG45A(56),SIG60A(56)
77*KEEP,NKGI.
78 COMMON /NKGI/ SEL,SELLG,STH,ZEL,ZELLG,ZSL,DIST,
79 * DISX,DISY,DISXY,DISYX,DLAX,DLAY,DLAXY,DLAYX,
80 * OBSATI,RADNKG,RMOL,TLEV,TLEVCM,IALT
81 DOUBLE PRECISION SEL(10),SELLG(10),STH(10),ZEL(10),ZELLG(10),
82 * ZSL(10),DIST(10),
83 * DISX(-10:10),DISY(-10:10),
84 * DISXY(-10:10,2),DISYX(-10:10,2),
85 * DLAX (-10:10,2),DLAY (-10:10,2),
86 * DLAXY(-10:10,2),DLAYX(-10:10,2),
87 * OBSATI(2),RADNKG,RMOL(2),TLEV(10),TLEVCM(10)
88 INTEGER IALT(2)
89*KEEP,OBSPAR.
90 COMMON /OBSPAR/ OBSLEV,THCKOB,XOFF,YOFF,THETAP,PHIP,
91 * THETPR,PHIPR,NOBSLV
92 DOUBLE PRECISION OBSLEV(10),THCKOB(10),XOFF(10),YOFF(10),
93 * THETAP,THETPR(2),PHIP,PHIPR(2)
94 INTEGER NOBSLV
95*KEEP,PAM.
96 COMMON /PAM/ PAMA,SIGNUM
97 DOUBLE PRECISION PAMA(6000),SIGNUM(6000)
98*KEEP,PARPAR.
99 COMMON /PARPAR/ CURPAR,SECPAR,PRMPAR,OUTPAR,C,
100 * E00,E00PN,PTOT0,PTOT0N,THICKH,ITYPE,LEVL
101 DOUBLE PRECISION CURPAR(14),SECPAR(14),PRMPAR(14),OUTPAR(14),
102 * C(50),E00,E00PN,PTOT0,PTOT0N,THICKH
103 INTEGER ITYPE,LEVL
104*KEEP,PARPAE.
105 DOUBLE PRECISION GAMMA,COSTHE,PHI,H,T,X,Y,CHI,BETA,GCM,ECM
106 EQUIVALENCE (CURPAR(2),GAMMA), (CURPAR(3),COSTHE),
107 * (CURPAR(4), PHI ), (CURPAR(5), H ),
108 * (CURPAR(6), T ), (CURPAR(7), X ),
109 * (CURPAR(8), Y ), (CURPAR(9), CHI ),
110 * (CURPAR(10),BETA), (CURPAR(11),GCM ),
111 * (CURPAR(12),ECM )
112*KEEP,PRIMSP.
113 COMMON /PRIMSP/ PSLOPE,LLIMIT,ULIMIT,LL,UL,SLEX,ISPEC
114 DOUBLE PRECISION PSLOPE,LLIMIT,ULIMIT,LL,UL,SLEX
115 INTEGER ISPEC
116*KEEP,RANDPA.
117 COMMON /RANDPA/ FAC,U1,U2,RD,NSEQ,ISEED,KNOR
118 DOUBLE PRECISION FAC,U1,U2
119 REAL RD(3000)
120 INTEGER ISEED(103,10),NSEQ
121 LOGICAL KNOR
122*KEEP,RANGE.
123 COMMON /RANGE/ CC
124 DOUBLE PRECISION CC(20)
125*KEEP,RECORD.
126 COMMON /RECORD/ IRECOR
127 INTEGER IRECOR
128*KEEP,RUNPAR.
129 COMMON /RUNPAR/ FIXHEI,THICK0,HILOECM,HILOELB,
130 * STEPFC,NRRUN,NSHOW,PATAPE,MONIIN,
131 * MONIOU,MDEBUG,NUCNUC,
132 * CETAPE,
133 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
134 * N1STTR,MDBASE,
135 * DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
136 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
137 * ,GHEISH,GHESIG
138 COMMON /RUNPAC/ DSN,HOST,USER
139 DOUBLE PRECISION FIXHEI,THICK0,HILOECM,HILOELB
140 REAL STEPFC
141 INTEGER NRRUN,NSHOW,PATAPE,MONIIN,MONIOU,MDEBUG,NUCNUC,
142 * SHOWNO,ISHW,NOPART,NRECS,NBLKS,MAXPRT,NDEBDL,
143 * N1STTR,MDBASE
144 INTEGER CETAPE
145 CHARACTER*79 DSN
146 CHARACTER*20 HOST,USER
147
148 LOGICAL DEBDEL,DEBUG,FDECAY,FEGS,FIRSTI,FIXINC,FIXTAR,
149 * FIX1I,FMUADD,FNKG,FPRINT,FDBASE
150 * ,GHEISH,GHESIG
151*KEEP,STACKF.
152 COMMON /STACKF/ STACK,STACKP,EXST,NSHIFT,NOUREC,ICOUNT,NTO,NFROM
153 INTEGER MAXSTK
154 PARAMETER (MAXSTK = 12*340*2)
155 DOUBLE PRECISION STACK(MAXSTK)
156 INTEGER STACKP,EXST,NSHIFT,NOUREC,ICOUNT,NTO,NFROM
157*KEEP,STRBAR.
158 COMMON /STRBAR/ CSTRBA
159 DOUBLE PRECISION CSTRBA(11)
160*KEEP,VERS.
161 COMMON /VERS/ VERNUM,MVDATE,VERDAT
162 DOUBLE PRECISION VERNUM
163 INTEGER MVDATE
164 CHARACTER*18 VERDAT
165*KEEP,VENUS.
166 COMMON /VENUS/ ISH0,IVERVN,MTAR99,FVENUS,FVENSG
167 INTEGER ISH0,IVERVN,MTAR99
168 LOGICAL FVENUS,FVENSG
169*KEEP,CEREN3.
170 COMMON /CEREN3/ CERCNT,DATAB2,LHCER
171 INTEGER MAXBF2
172 PARAMETER (MAXBF2 = 39 * 7)
173 DOUBLE PRECISION CERCNT
174 REAL DATAB2(MAXBF2)
175 INTEGER LHCER
176*KEND.
177
178 DOUBLE PRECISION COAN,SE,TEMP1,TEMP2,TEMP3,THICK,TTIME,ZE,ZS,ZX
179 INTEGER I,IA,J,L,N
180 EXTERNAL THICK
181 CHARACTER*1 MARK
182C-----------------------------------------------------------------------
183
184C SAY HELLO
185 WRITE(MONIOU,112)
186 112 FORMAT(/' ',120('A')//
187 *' OOO OOO OOOO OOOO OO O O O '/
188 *' O O O O O O O O OO O O O O '/
189 *' O O O O O O OO O O O O '/
190 *' O O O O O OOOO OO OO O O'/
191 *' O O O OOOO O OO O O OOOOOOO'/
192 *' O O O O O O O O OO O O O O'/
193 *' OOO OOO O O OOOO OO O O O O'//
194 *' COSMIC RAY SIMULATION FOR KASCADE'///
195 *' A PROGRAM TO SIMULATE EXTENSIVE AIR SHOWERS IN ATMOSPHERE'//
196 *' BASED ON A PROGRAM OF P.K.F. GRIEDER, UNIVERSITY BERN,',
197 *' SWITZERLAND'/
198 *' HDPM MODEL ACCORDING TO J.N. CAPDEVIELLE, COLLEGE DE FRANCE,',
199 *' PARIS, FRANCE'/
200 *' VENUS MODEL ACCORDING TO K. WERNER, UNIVERSITY NANTES, FRANCE'/
201 *' GHEISHA ROUTINES ACCORDING TO H. FESEFELDT, RWTH. AACHEN,',
202 *' GERMANY'/
203 *' EGS4 AND NKG FORMULAS FOR SIMULATION OF EL.MAG. PARTICLES'//)
204
205 MARK = '1'
206
207 WRITE(MONIOU,912) VERNUM,MARK,VERDAT
208 912 FORMAT(' INSTITUT FUER KERNPHYSIK '/
209 * ' FORSCHUNGSZENTRUM UND UNIVERSITAET KARLSRUHE'/
210 * ' POSTFACH 3640'/
211 * ' D-76021 KARLSRUHE'/
212 * ' GERMANY'//
213 * ' IN CASE OF PROBLEMS CONTACT:'/
214 * ' DIETER HECK JOHANNES KNAPP'/
215 * ' E-MAIL: HECK@IK3.FZK.DE KNAPP@IK1.FZK.DE'/
216 * ' FAX: (49) 7247-82-4075 (49) 7247-82-3548'/
217 * ' TEL: (49) 7247-82-3777 (49) 7247-82-3549'//
218 * ' NUMBER OF VERSION : ',F6.3,A1/
219 * ' DATE OF VERSION : ',A18 /)
220
221 WRITE(MONIOU,141)
222 141 FORMAT(//' CERENKOV RADIATION IS GENERATED'/
223 * ' ==============================='//)
224
225C INITIALIZE FIELD WITH PARTICLE MASSES
226 CALL PAMAF
227
228
229C READ RUN STEERING DATA CARDS
230 CALL DATAC
231
232C CLEARS BUFFERS FOR HEADER AND FILLS IN PERMANENT INFORMATION
233 DO 889 L = 1,MAXBUF
234 EVTH(L) = 0.
235 EVTE(L) = 0.
236 RUNH(L) = 0.
237 RUNE(L) = 0.
238 DATAB(L) = 0.
239 DATAB2(L) = 0.
240 889 CONTINUE
241
242
243C PERMANENT INFORMATION
244C CHARACTER STRINGS
245 CRUNH = 'RUNH'
246 CRUNE = 'RUNE'
247 CEVTH = 'EVTH'
248 CEVTE = 'EVTE'
249
250 RUNH(2) = NRRUN
251 RUNE(2) = NRRUN
252 EVTH(44) = NRRUN
253
254C DATE OF RUN
255 WRITE(MONIOU,101)
256 101 FORMAT(//' ',10('='),' START OF RUN ',55('='))
257 CALL PRTIME(TTIME)
258 RUNH(3) = TTIME
259 EVTH(45) = TTIME
260
261C VERSION OF PROGRAM
262 RUNH(4) = VERNUM
263 EVTH(46) = VERNUM
264
265C-----------------------------------------------------------------------
266C INITIALISATION FOR RANDOM NUMBER GENERATOR
267 IF ( FEGS .AND. NSEQ .LT. 2 ) NSEQ = 2
268C CERENKOV SELECTION DEMANDS ALWAYS EGS CALCULATION
269 FEGS = .TRUE.
270C IN CASE OF CERENKOV CALCULATIONS THE 3. RANDOM SEQUENCE IS NEEDED
271 IF ( NSEQ .LT. 3 ) NSEQ = 3
272 DO 281 I = 1,NSEQ
273 IF ( .NOT. DEBUG .AND. .NOT. DEBDEL .AND.
274 * (ISEED(2,I) .GT. 1000 .OR. ISEED(3,I) .GT. 0) ) THEN
275 WRITE(MONIOU,2811) I
2762811 FORMAT(/' #########################################'/
277 * ' ## IMPROPER INITIALIZATION OF RANDOM ##'/
278 * ' ## NUMBER GENERATOR SEQUENCE ',I6,' ##'/
279 * ' ## IS EXTREMELY TIME CONSUMING ##'/
280 * ' ## PLEASE READ THE MANUALS ##'/
281 * ' #########################################'/)
282 ENDIF
283 CALL RMMAQ( ISEED(1,I), I, 'S' )
284 281 CONTINUE
285 KNOR = .TRUE.
286
287 WRITE(MONIOU,158) (L,(ISEED(J,L),J=1,3),L=1,NSEQ)
288 158 FORMAT (/' RANDOM NUMBER GENERATOR AT BEGIN OF RUN :'/
289 * (' SEQUENCE = ',I2,' SEED = ',I9,' CALLS = ',I9,
290 * ' BILLIONS = ',I9))
291
292C-----------------------------------------------------------------------
293C READ CROSS SECTIONS AND PROBABILITIES FOR NUCLEUS-NUCLEUS COLLISIONS
294 OPEN(UNIT=NUCNUC,FILE='NUCNUCCS',STATUS='OLD')
295 READ(NUCNUC,500) SIGN30,SIGN45,SIGN60,SIGO30,SIGO45,SIGO60,
296 * SIGA30,SIGA45,SIGA60
297 READ(NUCNUC,500) (PNOA30(I,1),I=1,1540),(PNOA45(I,1),I=1,1540),
298 * (PNOA60(I,1),I=1,1540),(PNOA30(I,2),I=1,1540),
299 * (PNOA45(I,2),I=1,1540),(PNOA60(I,2),I=1,1540),
300 * (PNOA30(I,3),I=1,1540),(PNOA45(I,3),I=1,1540),
301 * (PNOA60(I,3),I=1,1540)
302 500 FORMAT( 5E16.10 )
303 CLOSE(UNIT=NUCNUC)
304
305C INELASTIC CROSS SECTIONS FOR PROJECTICLE WITH MASS NUMBER IA
306 DO 501 IA = 1,56
307 SIG30A(IA) = COMPOS(1)*SIGN30(IA) + COMPOS(2)*SIGO30(IA)
308 * + COMPOS(3)*SIGA30(IA)
309 SIG45A(IA) = COMPOS(1)*SIGN45(IA) + COMPOS(2)*SIGO45(IA)
310 * + COMPOS(3)*SIGA45(IA)
311 SIG60A(IA) = COMPOS(1)*SIGN60(IA) + COMPOS(2)*SIGO60(IA)
312 * + COMPOS(3)*SIGA60(IA)
313
314 IF (DEBUG) WRITE(MDEBUG,544) IA,SIG30A(IA),SIG45A(IA),SIG60A(IA)
315 544 FORMAT(' START : CROSS SECTIONS A-AIR : A=',I2,1P,3E14.6)
316 501 CONTINUE
317
318 WRITE(MONIOU,503)
319 503 FORMAT (//' ',10('='),' INTERACTION MODELS ',49('='))
320C HIGH ENERGY HADRONIC INTERACTION MODEL
321 IF ( FVENUS ) THEN
322 WRITE(MONIOU,*) 'VENUS TREATS HIGH ENERGY HADRONIC INTERACTIONS'
323 CALL VENINI
324 IF ( .NOT. GHEISH ) THEN
325 GHEISH = .TRUE.
326 WRITE(MONIOU,*)'GHEISHA OPTION NOT SELECTED, BUT SWITCHED ON'
327 ENDIF
328 IF ( NFRAGM .EQ. 0 ) THEN
329 WRITE(MONIOU,*)
330 * ' TOTAL FRAGMENTATION OF PRIMARY NUCLEUS IN FIRST INTERACTION'
331 ELSEIF ( NFRAGM .EQ. 1 ) THEN
332 WRITE(MONIOU,*)
333 * ' NO FRAGMENTATION, NO EVAPORATION OF REMAINDER'
334 ELSEIF ( NFRAGM .EQ. 2 ) THEN
335 WRITE(MONIOU,1504)
336 ELSEIF ( NFRAGM .EQ. 3 ) THEN
337 WRITE(MONIOU,1505)
338 ELSE
339 NFRAGM = 4
340 WRITE(MONIOU,1507)
341 ENDIF
342 WRITE(MONIOU,*)
343 ELSE
344 WRITE(MONIOU,1506)
345 ENDIF
3461506 FORMAT(' HDPM ROUTINES TREAT HIGH ENERGY HADRONIC INTERACTIONS')
347
348
349 IF ( .NOT. FVENUS ) THEN
350C INPUT FLAGS FOR HDPM OPTIONS
351 WRITE(MONIOU,*)'HDPM GENERATOR SPECIFICATIONS ARE:'
352 IF ( NFLAIN .EQ. 0 ) THEN
353 WRITE(MONIOU,*) ' RANDOM NUMBER OF INTERACTIONS IN AIR TARGET'
354 IF ( NFLDIF .EQ. 0 ) THEN
355 WRITE(MONIOU,*) ' NO DIFFRACTIVE SECOND INTERACTIONS'
356 ELSE
357 WRITE(MONIOU,*) ' DIFFRACTIVE SECOND INTERACTIONS'
358 ENDIF
359 ELSE
360 WRITE(MONIOU,*) ' FIXED NUMBER OF INTERACTIONS IN AIR TARGET'
361 ENDIF
362 IF ( NFLPI0 .EQ. 0 ) THEN
363 WRITE(MONIOU,*) ' RAPIDITY OF PI0 ACCORDING TO COLLIDER DATA'
364 ELSE
365 WRITE(MONIOU,*) ' RAPIDITY OF PI0 SAME AS THAT OF CHARGED'
366 ENDIF
367 IF ( NFLPIF .EQ. 0 ) THEN
368 WRITE(MONIOU,*) ' NO FLUCTUATIONS OF NUMBER OF PI0'
369 ELSE
370 WRITE(MONIOU,*)' FLUCTUATIONS OF NUMBER OF PI0 AS MEASURED ',
371 * 'AT THE COLLIDER'
372 ENDIF
373 IF ( NFLCHE .EQ. 0 ) THEN
374 WRITE(MONIOU,*) ' CHARGE EXCHANGE INTERACTION POSSIBLE '
375 ELSE
376 WRITE(MONIOU,*) ' NO CHARGE EXCHANGE INTERACTION POSSIBLE '
377 ENDIF
378 IF ( NFRAGM .EQ. 0 ) THEN
379 WRITE(MONIOU,*)' TOTAL FRAGMENTION OF PRIMARY NUCLEUS IN ',
380 * 'FIRST INTERACTION'
381 ELSEIF ( NFRAGM .EQ. 1 ) THEN
382 WRITE(MONIOU,*) ' NO FRAGMENTATION, NO EVAPORATION OF REMAINDER'
383 ELSEIF ( NFRAGM .EQ. 2 ) THEN
384 WRITE(MONIOU,1504)
3851504 FORMAT(' NO FRAGMENTATION, EVAPORATION OF REMAINDER ',
386 * ' (PT AFTER JACEE)')
387 ELSEIF ( NFRAGM .EQ. 3 ) THEN
388 WRITE(MONIOU,1505)
3891505 FORMAT(' NO FRAGMENTATION, EVAPORATION OF REMAINDER ',
390 * ' (PT AFTER GOLDHABER)')
391 ELSE
392 NFRAGM = 4
393 WRITE(MONIOU,1507)
3941507 FORMAT(' NO FRAGMENTATION, EVAPORATION OF REMAINDER ',
395 * ' (WITH PT = 0.)')
396 ENDIF
397 ENDIF
398 WRITE(MONIOU,*)
399
400C LOW ENERGY HADRONIC INTERACTION MODEL
401 IF ( GHEISH ) THEN
402 WRITE(MONIOU,*) 'GHEISHA TREATS LOW ENERGY HADRONIC ',
403 * 'INTERACTIONS'
404 CALL CGHINI
405 ELSE
406 WRITE(MONIOU,*) 'ISOBAR ROUTINES TREAT LOW ENERGY HADRONIC ',
407 * 'INTERACTIONS'
408 HILOELB = 53.D0
409 ENDIF
410
411C WRITE HADRONIC STEERING FLAGS TO RUNHEADER
412 RUNH(270) = NFLAIN
413 RUNH(271) = NFLDIF
414 RUNH(272) = NFLPI0 + 100. * NFLPIF
415 RUNH(273) = NFLCHE + 100. * NFRAGM
416
417 EVTH(65) = NFLAIN
418 EVTH(66) = NFLDIF
419 EVTH(67) = NFLPI0
420 EVTH(68) = NFLPIF
421 EVTH(69) = NFLCHE
422 EVTH(70) = NFRAGM
423
424 HILOECM = SQRT(2.D0*PAMA(14)*(PAMA(14) + HILOELB))
425 WRITE(MONIOU,*) 'START: HIGH ENERGY INTERACTION MODEL USED ABOVE'
426 WRITE(MONIOU,*) ' ',HILOELB,' GEV LAB ENERGY OR'
427 WRITE(MONIOU,*) ' ',HILOECM,' GEV CM ENERGY'
428
429C INPUT STEERING FLAGS FOR ELECTROMAGNETIC PART
430 WRITE(MONIOU,*)
431 IF ( FNKG ) THEN
432 WRITE(MONIOU,*)'ELECTROMAGNETIC COMPONENT SIMULATED WITH NKG'
433 IF ( ULIMIT .GT. 2.D7 ) THEN
434 WRITE(MONIOU,*)'#############################################'
435 WRITE(MONIOU,*)'# W A R N I N G NKG IS WITHOUT LPM EFFECT #'
436 WRITE(MONIOU,*)'#############################################'
437 ENDIF
438 WRITE(MONIOU,*)
439 ENDIF
440 IF ( FEGS ) THEN
441 WRITE(MONIOU,*)'ELECTROMAGNETIC COMPONENT SIMULATED WITH EGS4'
442 WRITE(MONIOU,*)
443 ENDIF
444 IF ( .NOT. (FNKG .OR. FEGS) ) WRITE(MONIOU,*)
445 * 'ELECTROMAGNETIC COMPONENT IS NOT SIMULATED'
446 IF ( FEGS ) THEN
447 IF ( STEPFC .GT. 10. .OR. STEPFC .LE. 0. ) THEN
448 WRITE(MONIOU,*)'STEP LENGTH FACTOR FOR ELECTRON MULTIPLE ',
449 * 'SCATTERING =',STEPFC,' NOT CORRECT'
450 WRITE(MONIOU,*)'PLEASE READ THE MANUALS'
451 STOP
452 ENDIF
453 IF ( STEPFC .LT. 10. ) WRITE(MONIOU,*)'STEP LENGTH ',
454 * 'FACTOR FOR ELECTRON MULTIPLE SCATTERING =',STEPFC
455C INITIALIZE EGS4 PACKAGE
456 CALL EGSINI
457 IF ( ULIMIT .GT. 2.D7 ) THEN
458 WRITE(MONIOU,*)'#############################################'
459 WRITE(MONIOU,*)'# W A R N I N G EGS IS WITHOUT LPM EFFECT #'
460 WRITE(MONIOU,*)'#############################################'
461 ENDIF
462 ENDIF
463C WRITE STEERING FLAGS FOR ELECTROMAGNETIC PART AS REAL TO HEADER
464 IF ( FNKG ) THEN
465 RUNH(20) = 1.
466 EVTH(74) = 1.
467 ELSE
468 RUNH(20) = 0.
469 EVTH(74) = 0.
470 ENDIF
471 IF ( FEGS ) THEN
472 RUNH(19) = 1.
473 EVTH(73) = 1.
474 ELSE
475 RUNH(19) = 0.
476 EVTH(73) = 0.
477 ENDIF
478
479 EVTH(95) = STEPFC
480
481C PROGRAM CONFIGURATIONS FOR EVENT HEADER
482 IF ( GHEISH ) THEN
483 EVTH(75) = 1.
484 ELSE
485 EVTH(75) = 0.
486 ENDIF
487 IF ( FVENUS ) THEN
488 EVTH(76) = 1.
489 ELSE
490 EVTH(76) = 0.
491 ENDIF
492 EVTH(139) = 0.
493 EVTH(140) = 0.
494 EVTH(141) = 0.
495 EVTH(142) = 0.
496 EVTH(143) = 0.
497 EVTH(144) = 0.
498 EVTH(145) = 0.
499 EVTH(77) = 1.
500 EVTH(78) = 0.
501 EVTH(79) = 0.
502 EVTH(80) = 3.
503
504C-----------------------------------------------------------------------
505C BEGIN OF TAPE FOR IBM, FOR TRANSPUTER SEE BEGIN OF EVT
506
507C-----------------------------------------------------------------------
508C PHYSICAL CONSTANTS
509 ENEPER = EXP(1.D0)
510 C(6) = ( PAMA(5) / PAMA(11) )**2
511 C(7) = ( PAMA(5) / PAMA(8) )**2
512 C(8) = ( PAMA(5)**2 + PAMA(2)**2 ) * 0.5D0 / PAMA(5)
513 C(20) = 10.D0 * C(21)
514 C(27) = COS( C(26) )
515 C(29) = COS( C(28) )
516 C(44) = MAX( PAMA(8)+C(4), PAMA(14)+C(5) )
517 C(45) = PAMA(8) * PAMA(14) * 2.D0
518 C(46) = PAMA(8)**2 + PAMA(14)**2
519 C(48) = (PAMA(8)**2 + PAMA(5)**2) / (2.D0*PAMA(8)*PAMA(5))
520 C(49) = SQRT(C(48)**2 - 1.D0) / C(48)
521
522 CKA(13) = 2.D0 * PAMA(11) * PAMA(14)
523 CKA(14) = PAMA(11)**2 + PAMA(14)**2
524 CKA(17) = SQRT( ( (PAMA(11)**2 + PAMA(5)**2)
525 * / (2.D0*PAMA(11)) )**2 - PAMA(5)**2 )
526 CKA(18) = SQRT( ( (PAMA(11)**2 + PAMA(8)**2 - PAMA(7)**2)
527 * / (2.D0*PAMA(11)) )**2 - PAMA(8)**2 )
528 CKA(22) = MAX( C(5)+PAMA(14), PAMA(11)+C(4) )
529 CKA(28) = SQRT(1.D0 + CKA(17)**2/PAMA(5)**2)
530 CKA(29) = SQRT(1.D0 - 1.D0/CKA(28)**2)
531 CKA(30) = SQRT(1.D0 + CKA(18)**2/PAMA(8)**2)
532 CKA(31) = SQRT(1.D0 - 1.D0/CKA(30)**2)
533 CKA(41) = PAMA(16)
534 CKA(42) = (PAMA(11)**2 + PAMA(7)**2 - PAMA(8)**2) /
535 * (2.D0*PAMA(11)*PAMA(7))
536 CKA(43) = CKA(41) / (2.D0*PAMA(7))
537 CKA(44) = SQRT(1.D0 - 1.D0/CKA(43)**2)
538 CKA(45) = CKA(41) / (2.D0*PAMA(8))
539 CKA(46) = SQRT(1.D0 - 1.D0/CKA(45)**2)
540
541C SET CONSTANTS FOR MUON BREMSSTRAHLUNG
542 CMUON(3) = 7.D0**OB3
543 CMUON(6) = 8.D0**OB3
544 CMUON(9) = 18.D0**OB3
545 CMUON(1) = LOG( 189.D0 * PAMA(5) / (CMUON(3)*PAMA(2)) )
546 CMUON(4) = LOG( 189.D0 * PAMA(5) / (CMUON(6)*PAMA(2)) )
547 CMUON(7) = LOG( 189.D0 * PAMA(5) / (CMUON(9)*PAMA(2)) )
548 * + LOG( TB3/CMUON(9) )
549 SE = SQRT(EXP(1.D0))
550 CMUON(2) = 189.D0 * SE*PAMA(5)**2/(2.D0*PAMA(2)*CMUON(3))
551 CMUON(5) = 189.D0 * SE*PAMA(5)**2/(2.D0*PAMA(2)*CMUON(6))
552 CMUON(8) = 189.D0 * SE*PAMA(5)**2/(2.D0*PAMA(2)*CMUON(9))
553 CMUON(10) = 0.75D0 * PAMA(5) * SE
554 CMUON(3) = CMUON(3) * CMUON(10)
555 CMUON(6) = CMUON(6) * CMUON(10)
556 CMUON(9) = CMUON(9) * CMUON(10)
557 CMUON(11) = LOG( BCUT/PAMA(5) )
558
559 DO 1 I = 1,50
560 CANN(I) = 0.D0
561 1 CONTINUE
562 COAN = 0.D0
563 DO 25 N = 1,12
564 COAN = COAN + CAN(N)
565 CANN(N) = COAN
566 25 CONTINUE
567 COAN = 0.D0
568 DO 26 N = 13,26
569 COAN = COAN + CAN(N)
570 CANN(N) = COAN
571 26 CONTINUE
572
573C-----------------------------------------------------------------------
574C INITIALIZE CONSTANTS FOR MUON MULTIPLE SCATTERING (MOLIERE)
575C SEE SUBROUTINE GMOLI OF GEANT321 (CERN)
576 IF (FMOLI) THEN
577 TEMP1 = COMPOS(1) * 7.D0 * 8.D0 / 14.D0
578 TEMP2 = COMPOS(2) * 8.D0 * 9.D0 / 16.D0
579 TEMP3 = COMPOS(3) * 18.D0 * 19.D0 / 40.D0
580 ZS = TEMP1 + TEMP2 + TEMP3
581 ZE = -TB3*(TEMP1*LOG(7.D0) +TEMP2*LOG(8.D0) +TEMP3*LOG(18.D0))
582 ZX = TEMP1*LOG(1.D0 + 3.34D0 * ( 7.D0/C(50))**2)
583 * +TEMP2*LOG(1.D0 + 3.34D0 * ( 8.D0/C(50))**2)
584 * +TEMP3*LOG(1.D0 + 3.34D0 * (18.D0/C(50))**2)
585C NOTE: CHC IS DEFINED DIFFERENT FROM GEANT WITHOUT DENSITY
586 CHC = 0.39612D-3 * SQRT(ZS)
587C NOTE: OMC IS DEFINED DIFFERENT FROM GEANT WITHOUT DENSITY
588 OMC = 6702.33D0 * ZS * EXP( (ZE-ZX)/ZS )
589 EVTH(146) = 1.
590 ELSE
591 EVTH(146) = 0.
592 ENDIF
593
594C-----------------------------------------------------------------------
595C TEST ON INPUT VALUES
596
597C PRINT CONTROL OUTPUT
598 IF ( CC(1) .GE. CC(2) .OR.
599 * CC(2) .GE. CC(3) .OR.
600 * CC(3) .GE. CC(4) .OR.
601 * CC(5) .GE. CC(6) .OR.
602 * CC(6) .GE. CC(7) .OR.
603 * CC(7) .GE. CC(8) .OR.
604 * CC(9) .GE. CC(10) .OR.
605 * CC(10) .GE. CC(11) .OR.
606 * CC(11) .GE. CC(12) .OR.
607 * PAMA(14)+C(3) .GT. CC(1) .OR.
608 * PAMA(14)+C(4) .GT. CC(2) .OR.
609 * C(4)*2. .GT. CC(3) .OR.
610 * C(3)+PAMA(8) .GT. CC(5) .OR.
611 * C(44) .GT. CC(6) .OR.
612 * C(4)+C(5) .GT. CC(7) .OR.
613 * PAMA(14)+C(4) .GE. C(4)*2. .OR.
614 * C(44) .GE. C(4)+C(5) ) THEN
615 WRITE(MONIOU,106)
616 106 FORMAT (' ERROR OR INCOMPATIBILITY IN CONSTANTS')
617C PRINT CONTROL OUTPUT
618 WRITE(MONIOU,103) (C(I),I=1,50)
619 103 FORMAT (//' ',10('='),' CONSTANTS AND PARAMETERS ',43('=')
620 * //' PHYSICAL CONSTANTS (C)' // (1P,4(E15.8,1X),E15.8) )
621 WRITE(MONIOU,110) (CKA(I),I=1,80)
622 110 FORMAT (//' CONSTANTS FOR KAONS CKA(1) TO CKA(40)'
623 * // (1P,4(E15.8,1X),E15.8) )
624 WRITE(MONIOU,114) (CETA(I),I=1,5)
625 114 FORMAT (//' CONSTANTS FOR ETAS CETA(1) TO CETA(5)'
626 * // (1P,4(E15.8,1X),E15.8) )
627 WRITE(MONIOU,115) (CSTRBA(I),I=1,11)
628 115 FORMAT (//' CONSTANTS FOR STRANGE BARYONS CSTRBA(1) TO ',
629 * 'CSTRBA(11)'// (1P,4(E15.8,1X),E15.8) )
630 WRITE(MONIOU,206) (CAN(I),I=1,30)
631 206 FORMAT (//' ANNIHILATION PARAMETERS, SET 1 (CAN)'
632 * // (1P,4(E15.8,1X),E15.8) )
633 WRITE(MONIOU,209) (CANN(I),I=1,30)
634 209 FORMAT (//' ANNIHILATION PARAMETERS, SET 2 (CANN)'
635 * // (1P,4(E15.8,1X),E15.8) )
636 WRITE(MONIOU,60) (CC(I),I=1,12)
637 60 FORMAT (//' THRESHOLD ENERGIES OF INTERACTION INTERVALS IN GEV',
638 * ' (CC)'// (1P,4(E15.8,1X),E15.8) )
639 WRITE(MONIOU,106)
640 STOP
641 ENDIF
642
643C FILL CONSTANTS IN RUN HEADER
644 DO 3001 L = 1,50
645 RUNH(24+L) = C(L)
646 RUNH(154+L) = CAN(L)
647 RUNH(204+L) = CANN(L)
648 3001 CONTINUE
649 DO 3002 L = 1,20
650 RUNH(74+L) = CC(L)
651 3002 CONTINUE
652 DO 3003 L = 1,40
653 RUNH(94+L) = CKA(L)
654 3003 CONTINUE
655 DO 3004 L = 1,5
656 RUNH(134+L) = CETA(L)
657 3004 CONTINUE
658 DO 3005 L = 1,11
659 RUNH(139+L) = CSTRBA(L)
660 3005 CONTINUE
661 DO 3007 L = 1,5
662 RUNH(254+L) = AATM(L)
663 RUNH(259+L) = BATM(L)
664 RUNH(264+L) = CATM(L)
665 DATM(L) = 1.D0 / CATM(L)
666 3007 CONTINUE
667
668C SET LOWER BOUNDARIES OF THE AIR LAYERS
669 HLAY(1) = 0.D0
670 HLAY(2) = 4.D5
671 HLAY(3) = 1.D6
672 HLAY(4) = 4.D6
673 HLAY(5) = 1.D7
674C CALCULATE THICKNESS AT LOWER BOUNDARIES OF AIR LAYERS
675 DO 100 L= 1,5
676 THICKL(L) = THICK(HLAY(L))
677 100 CONTINUE
678
679 CALL STAEND
680
681 RETURN
682 END
Note: See TracBrowser for help on using the repository browser.