fe4da5cc |
1 | *CMZ : 19/11/98 21.32.03 by Federico Carminati |
2 | *-- Author : |
3 | C*********************************************************************** |
4 | C |
5 | SUBROUTINE ZDC_ADDANG(THETA1,PHI1,THETA2,PHI2,THETA3,PHI3) |
6 | C |
7 | C This subroutine takes spherical polar angles THETA1 and PHI1 |
8 | C and adds to them the direction PHI2 on a cone of half-opening angle |
9 | C THETA2 (= angle from cone axis to cone surface) to produce |
10 | C the new space angles THETA3 and PHI3. All angles are in |
11 | C degrees, with 0<=THETA<=180 and 0<=PHI<360. There is no |
12 | C other restriction on THETA1, PHI1, THETA2 or PHI2. |
13 | C |
14 | C 27-AUG-1992 S. Coutu |
15 | C |
16 | C*********************************************************************** |
17 | #undef CERNLIB_GEANT321_GCONSP_INC |
18 | #include "geant321/gconsp.inc" |
19 | DOUBLE PRECISION THETA1,PHI1,THETA2,PHI2,THETA3,PHI3 |
20 | DOUBLE PRECISION TEMP,CX,CY,CZ,CT1,ST1,CT2,ST2, |
21 | + CP1,SP1,CP2,SP2,RTHETA3 |
22 | C |
23 | CT1=COS(THETA1*DEGRAD) |
24 | ST1=SIN(THETA1*DEGRAD) |
25 | CP1=COS(PHI1*DEGRAD) |
26 | SP1=SIN(PHI1*DEGRAD) |
27 | CT2=COS(THETA2*DEGRAD) |
28 | ST2=SIN(THETA2*DEGRAD) |
29 | CP2=COS(PHI2*DEGRAD) |
30 | SP2=SIN(PHI2*DEGRAD) |
31 | CX=CT1*CP1*ST2*CP2+ST1*CP1*CT2-SP1*ST2*SP2 |
32 | CY=CT1*SP1*ST2*CP2+ST1*SP1*CT2+CP1*ST2*SP2 |
33 | CZ=CT1*CT2-ST1*ST2*CP2 |
34 | C |
35 | RTHETA3=ACOS(CZ) |
36 | THETA3=RTHETA3*RADDEG |
37 | IF(THETA3.EQ.0.OR.THETA3.EQ.180.)THEN |
38 | PHI3=0. |
39 | ELSE |
40 | TEMP=CX/SIN(RTHETA3) |
41 | IF(TEMP.GT.1.)TEMP=1. |
42 | IF(TEMP.LT.-1.)TEMP=-1. |
43 | PHI3=ACOS(TEMP)*RADDEG |
44 | IF(CY.LT.0.)PHI3=360.-PHI3 |
45 | ENDIF |
46 | * |
47 | END |
48 | |
49 | *CMZ : 21/11/98 17.33.14 by Federico Carminati |
50 | *-- Author : |
51 | C================================================================== |
52 | SUBROUTINE ZDC_EFERMI(ID) |
53 | * |
54 | #include "zdc_common.inc" |
55 | #undef CERNLIB_GEANT321_GCONSP_INC |
56 | #include "geant321/gconsp.inc" |
57 | DIMENSION RNDM(3) |
58 | * |
59 | CALL GRNDM(RNDM,3) |
60 | * |
61 | IF(ID.EQ.14)THEN |
62 | DO I=1,200 |
63 | IF(RNDM(3).GE.PROBINTP(I-1).AND.RNDM(3).LT.PROBINTP(I))GO TO |
64 | $ 10 |
65 | END DO |
66 | ELSE IF(ID.EQ.13) THEN |
67 | DO I=1,200 |
68 | IF(RNDM(3).GE.PROBINTN(I-1).AND.RNDM(3).LT.PROBINTN(I))GO TO |
69 | $ 10 |
70 | END DO |
71 | ENDIF |
72 | * |
73 | 10 PIPPO=PP(I)+0.001 |
74 | PHI=TWOPI*RNDM(1) |
75 | COST=1.-2*RNDM(2) |
76 | TET=ACOS(COST) |
77 | DDP(1)=PIPPO*SIN(TET)*COS(PHI) |
78 | DDP(2)=PIPPO*SIN(TET)*SIN(PHI) |
79 | DDP(3)=PIPPO*COST |
80 | * |
81 | END |
82 | |
83 | *CMZ : 23/02/99 07.53.48 by Federico Carminati |
84 | *CMZ : 2.03/01 06/08/98 10.13.59 by Federico Carminati |
85 | *-- Author : E.SCOMPARIN, 13/05/1996. |
86 | SUBROUTINE ZDC_END |
87 | C |
88 | C *** TERMINATION OF THE ZDC SIMULATION AT THE END OF A RUN *** |
89 | C |
90 | C CALLED BY : SXEND |
91 | C ORIGIN : E. SCOMPARIN |
92 | C |
93 | C |
94 | * CALL HROUT(0,ICYCLE,' ') |
95 | * CALL HREND('ALIZ') |
96 | * CLOSE(1) |
97 | |
98 | PRINT 10000 |
99 | 10000 FORMAT(1H /1H ,36('*'),' ZDC_END ',36('*')/ |
100 | $ 1H ,20X,'ZDC simulation package terminated'/ |
101 | $ 1H ,80('*')) |
102 | * |
103 | END |
104 | |
105 | *CMZ : 23/02/99 07.53.48 by Federico Carminati |
106 | *CMZ : 2.03/01 07/08/98 09.48.28 by Federico Carminati |
107 | *-- Author : E.SCOMPARIN, 13/05/1996. |
108 | SUBROUTINE ZDC_EVE |
109 | C |
110 | C *** TERMINATION OF THE ZDC SIMULATION AFTER EACH EVENT *** |
111 | C |
112 | C CALLED BY : SXEVE |
113 | C ORIGIN : E. SCOMPARIN |
114 | C |
115 | C -- RETRIEVE HITS INFORMATION FOR THE CURRENT EVENT |
116 | #undef CERNLIB_GEANT321_GCFLAG_INC |
117 | #include "geant321/gcflag.inc" |
118 | #include "zdc_common.inc" |
119 | * |
120 | INTEGER NUMVSP(NVZP),NUMVSN(NVZN), |
121 | + NUMBP(NVZP,1000),NUMBN(NVZN,1000), |
122 | + ITRP(1000),ITRN(1000) |
123 | REAL HITP(NHZ,1000),HITN(NHZ,1000) |
124 | * |
125 | REAL ENEH(1000),XZPH(1000),YZPH(1000),EZPH(1000,NZPTX,NZPTY), |
126 | + XZNH(1000),YZNH(1000),EZNH(1000,NZNTX,NZNTY) |
127 | INTEGER NPAH(1000),SFLH(1000),NPPH(1000,NZPTX,NZPTY), |
128 | + NPNH(1000,NZNTX,NZNTY) |
129 | * |
130 | DATA NUMVSP/0,0/ |
131 | DATA NUMVSN/0,0/ |
132 | * |
133 | DO JZERO=1,1000 |
134 | DO KZERO=1,NHZ |
135 | HITP(KZERO,JZERO)=0 |
136 | HITN(KZERO,JZERO)=0 |
137 | ENDDO |
138 | ENDDO |
139 | * |
140 | DO I=1,1000 |
141 | NPAH(I)=0 |
142 | XZPH(I)=99. |
143 | YZPH(I)=99. |
144 | XZNH(I)=99. |
145 | YZNH(I)=99. |
146 | ENDDO |
147 | * |
148 | ITOLD=0 |
149 | * CALL GTREVE |
150 | * |
151 | C -- RETRIEVE HITS INFORMATION FOR THE CURRENT EVENT |
152 | IF(ISECF.EQ.1.OR.(ISECF.EQ.0.AND.MOD(IEVENT,1).EQ.0))THEN |
153 | WRITE(6,*)'EVENT NO. ',IEVENT |
154 | CALL GPHITS('*','*') |
155 | WRITE(6,*) |
156 | ENDIF |
157 | * |
158 | NHMAXP=1000 |
159 | NHMAXN=1000 |
160 | * |
161 | CALL GFHITS(ZSET,ZP1,NVZP,NHZ,NHMAXP,0, |
162 | + NUMVSP,ITRP,NUMBP,HITP,NHITP) |
163 | CALL GFHITS(ZSET,ZN1,NVZN,NHZ,NHMAXN,0, |
164 | + NUMVSN,ITRN,NUMBN,HITN,NHITN) |
165 | * |
166 | * Fill n-tuple |
167 | * |
168 | ITCTP=0 |
169 | DO I=1,NHITP |
170 | IF(ITRP(I).NE.ITCTP)THEN |
171 | ITCTP=ITRP(I) |
172 | NPAH(ITCTP)=HITP(1,I) |
173 | ENEH(ITCTP)=HITP(2,I) |
174 | XZPH(ITCTP)=HITP(3,I) |
175 | YZPH(ITCTP)=HITP(4,I) |
176 | SFLH(ITCTP)=INT(HITP(5,I)) |
177 | ENDIF |
178 | J1=NUMBP(1,I) |
179 | J2=NUMBP(2,I) |
180 | EZPH(ITCTP,J1,J2)=HITP(7,I) |
181 | NPPH(ITCTP,J1,J2)=HITP(6,I) |
182 | ENDDO |
183 | ITCTN=0 |
184 | DO I=1,NHITN |
185 | IF(ITRN(I).NE.ITCTN)THEN |
186 | ITCTN=ITRN(I) |
187 | NPAH(ITCTN)=HITN(1,I) |
188 | ENEH(ITCTN)=HITN(2,I) |
189 | XZNH(ITCTN)=HITN(3,I) |
190 | YZNH(ITCTN)=HITN(4,I) |
191 | SFLH(ITCTN)=INT(HITN(5,I)) |
192 | ENDIF |
193 | J1=NUMBN(1,I) |
194 | J2=NUMBN(2,I) |
195 | EZNH(ITCTN,J1,J2)=HITN(7,I) |
196 | NPNH(ITCTN,J1,J2)=HITN(6,I) |
197 | ENDDO |
198 | * |
199 | DO I=1,1000 |
200 | IF(NPAH(I).NE.0)THEN |
201 | NPA=NPAH(I) |
202 | ENE=ENEH(I) |
203 | XZP=XZPH(I) |
204 | YZP=YZPH(I) |
205 | SFL=SFLH(I) |
206 | DO J=1,NZPTX |
207 | DO K=1,NZPTY |
208 | EZP(J,K)=EZPH(I,J,K) |
209 | NPP(J,K)=NPPH(I,J,K) |
210 | ENDDO |
211 | ENDDO |
212 | XZN=XZNH(I) |
213 | YZN=YZNH(I) |
214 | DO J=1,NZNTX |
215 | DO K=1,NZNTY |
216 | EZN(J,K)=EZNH(I,J,K) |
217 | NPN(J,K)=NPNH(I,J,K) |
218 | ENDDO |
219 | ENDDO |
220 | * CALL HFNT(10) |
221 | ENDIF |
222 | ENDDO |
223 | * |
224 | END |
225 | |
226 | *CMZ : 21/11/98 17.33.14 by Federico Carminati |
227 | *-- Author : |
228 | C================================================================== |
229 | SUBROUTINE ZDC_GFERMI(A,Z) |
230 | * |
231 | #include "zdc_common.inc" |
232 | #undef CERNLIB_GEANT321_GCONSP_INC |
233 | #include "geant321/gconsp.inc" |
234 | * CALCOLO DELLA DISTRIBUZIONE DEI MOMENTI SECONDO |
235 | * LA DISTRIBUZIONE DOPPIA GAUSSIANA (ILINOV), UGUALE |
236 | * PER NEUTRONI E PROTONI |
237 | * |
238 | PROBINTP(0)=0.0 |
239 | PROBINTN(0)=0.0 |
240 | SIG1=0.113 |
241 | SIG2=0.250 |
242 | ALFA=0.18*((A/12.)**0.333) |
243 | XK=(4.*PI)/((1.+ALFA)*TWOPI**1.5) |
244 | * |
245 | DO I=1,200 |
246 | P=FLOAT(I)*0.005 |
247 | PP(I)=P |
248 | E1=(P*P)/(2.*SIG1*SIG1) |
249 | E2=(P*P)/(2.*SIG2*SIG2) |
250 | F1=EXP(-(E1)) |
251 | F2=EXP(-(E2)) |
252 | PROBP=XK*P*P*(F1/(SIG1**3)+F2/(SIG2**3))*0.005 |
253 | PROBINTP(I)=PROBINTP(I-1)+PROBP |
254 | PROBINTN(I)=PROBINTP(I) |
255 | END DO |
256 | * |
257 | END |
258 | |
259 | *CMZ : 23/02/99 07.53.48 by Federico Carminati |
260 | *CMZ : 2.03/01 25/08/98 23.39.54 by Federico Carminati |
261 | *-- Author : |
262 | C-- Author : E.SCOMPARIN, 07/05/1996 |
263 | SUBROUTINE ZDC_INIT |
264 | C |
265 | C *** INITIALISATION FOR THE ZDC SIMULATION *** |
266 | C |
267 | C CALLED BY : SXINIT |
268 | C ORIGIN : E. SCOMPARIN |
269 | C |
270 | C |
271 | #undef CERNLIB_GEANT321_GCFLAG_INC |
272 | #include "geant321/gcflag.inc" |
273 | #undef CERNLIB_GEANT321_GCKINE_INC |
274 | #include "geant321/gckine.inc" |
275 | C |
276 | #include "zdc_common.inc" |
277 | * |
278 | * Initialize COMMON block ZDC_CGEOM |
279 | * |
280 | DATA HDZN/4.0,4.0,50.0/ |
281 | DATA HDZP/8.0,8.0,75.0/ |
282 | C Coordinates of the center of the ZDC front face in the MRS |
283 | DATA ZNPOS/-0.5,0.,11613./ |
284 | DATA ZPPOS/-19.0,0.,11563./ |
285 | DATA FIZN/0.,0.01825,50.0/ |
286 | DATA FIZP/0.,0.01825,75.0/ |
287 | DATA GRZN/0.025,0.025,50.0/ |
288 | DATA GRZP/0.040,0.040,75.0/ |
289 | DATA NCEN/11,11,0/ |
290 | DATA NCEP/10,10,0/ |
291 | * |
292 | * Initialize COMMON block ZDC_HITS |
293 | * |
294 | * Data for set ZDC |
295 | DATA ZSET,ZN1,ZP1 /'ZDC ', 'ZN1 ','ZP1 '/ |
296 | DATA CHNMSZN/'ZNTX','ZN1 '/ |
297 | DATA CHNMSZP/'ZPTX','ZP1 '/ |
298 | DATA NBITSZN /2,2/ |
299 | DATA NBITSZP /2,2/ |
300 | DATA IDTYPZN,IDTYPZP /1, 2/ |
301 | DATA NWHI,NWDI /1000, 1000/ |
302 | * |
303 | * Data for Detectors ZN1,ZP1 |
304 | * |
305 | DATA CHNAMH /'IPAR','EPAR','XPAR','YPAR','SFLG','FELE','ENER'/ |
306 | DATA NBITSH /8,32,32,32,1,32,32/ |
307 | DATA ORIG /0.,0.,10.,10.,0.,0.,0./ |
308 | DATA FACT /1.,1.,10000.,10000.,1.,10000.,10000./ |
309 | DATA NHSUM /2/ |
310 | * |
311 | * Read tables of Cerenkov light produced in the fibers |
312 | * |
313 | OPEN(67,FILE='ZDC/light22620362207s', |
314 | + FORM='FORMATTED', STATUS='OLD') |
315 | OPEN(68,FILE='ZDC/light22620362208s', |
316 | + FORM='FORMATTED', STATUS='OLD') |
317 | OPEN(69,FILE='ZDC/light22620362209s', |
318 | + FORM='FORMATTED', STATUS='OLD') |
319 | OPEN(70,FILE='ZDC/light22620362210s', |
320 | + FORM='FORMATTED', STATUS='OLD') |
321 | DO 10 K=1,NALFA |
322 | READ (67,10000)(TABLE(1,K,J),J=1,NBE) |
323 | READ (68,10000)(TABLE(2,K,J),J=1,NBE) |
324 | READ (69,10000)(TABLE(3,K,J),J=1,NBE) |
325 | READ (70,10000)(TABLE(4,K,J),J=1,NBE) |
326 | 10 CONTINUE |
327 | 10000 FORMAT(7(1X,F10.5)) |
328 | CLOSE (67) |
329 | CLOSE (68) |
330 | CLOSE (69) |
331 | CLOSE (70) |
332 | |
333 | IFERMI=1 |
334 | * |
335 | * Set Fermi flag according to SWITCH 3 |
336 | IF(ISWIT(3).EQ.1) IFERMI=1 |
337 | * define the LHC energy/nucleon |
338 | EPERN=2760. |
339 | * |
340 | * Open data file to read in particles from event generators |
341 | IF(IHIJ.GE.1)THEN |
342 | OPEN(27,FILE=FILEH,STATUS='OLD') |
343 | ELSE IF(IVNU.GE.1)THEN |
344 | OPEN(27,FILE=FILEV,STATUS='OLD') |
345 | ENDIF |
346 | C GET FERMI MOMENTUM DISTRIBUTIONS |
347 | CALL ZDC_GFERMI(207.,82.) |
348 | C |
349 | * CALL HROPEN (1,'ALIZ',FILEA,'N',1024,ISTAT) |
350 | * |
351 | * CALL HBNT(10,'ALIZDC','D') |
352 | * CALL HBNAME(10,'ALI',NPA, |
353 | * +'NPA[0,112]:U,ENE,SFL[0,1]:U,'// |
354 | * +'XZP:R*4:16:[-99,99],YZP:R*4:16:[-99,99],'// |
355 | * +'EZP(2,2),NPP(2,2),'// |
356 | * +'XZN:R*4:16:[-99,99],YZN:R*4:16:[-99,99],'// |
357 | * +'EZN(2,2),NPN(2,2)') |
358 | * |
359 | PRINT 10100 |
360 | 10100 FORMAT(1H /1H ,36('*'),' ZDC_INIT ',36('*')/ |
361 | + 1H ,20X,'ZDC simulation package initialised'/ |
362 | + 1H ,80('*')) |
363 | |
364 | END |
365 | |
366 | *CMZ : 21/11/98 17.33.14 by Federico Carminati |
367 | *-- Author : |
368 | SUBROUTINE ZDC_KINE(NT) |
369 | * Generation of event kinematics * |
370 | * ----------------------------------------------------------------* |
371 | #undef CERNLIB_GEANT321_GCKINE_INC |
372 | #include "geant321/gckine.inc" |
373 | #undef CERNLIB_GEANT321_GCFLAG_INC |
374 | #include "geant321/gcflag.inc" |
375 | #include "zdc_common.inc" |
376 | * |
377 | DIMENSION KATT(10000,2),PATT(10000,4) |
378 | * |
379 | INTEGER ITA(49),ITB(49) |
380 | REAL PLAB(3),UBUF(1),UB(1) |
381 | DIMENSION VERTEX(3) |
382 | REAL BALP(4),DDDP(4),PLABS(4) |
383 | * |
384 | DATA UBUF/0./ |
385 | DATA ITA/22,-11,11,0,-13,13,111,211,-211,130,321,-321, |
386 | +2112,2212,-2212,310,221,3122,3222,3212,3112,3322,3312, |
387 | +3334,-2112,-3122,-3112,-3212,-3222,-3322,-3312,-3334, |
388 | +-15,15,411,-411,421,-421,431,-431,4122,0,0,0,0,0,0,0,0/ |
389 | DATA ITB/10,-12,12,0,-14,14,110,120,-120,230,130,-130, |
390 | +1220,1120,-1120,0,220,2130,1130,1230,2230,1330,2330,3331, |
391 | +-1220,-2130,-1130,-1230,-2230,-1330,-2330,-3331,-16,16, |
392 | +-240,240,-140,140,0,0,2140,0,0,0,0,0,0,0,0/ |
393 | DATA IVR/0/ |
394 | * |
395 | DIMENSION RNDM(2) |
396 | CHARACTER*20 PARTNAME |
397 | * |
398 | * --- BEAMTY=1 Gaussian beam (width sigx sigy, mean values fx fy) |
399 | * --- BEAMTY=2 Uniform beam (width sigx sigy, mean values fx fy) |
400 | * Data cards BMTY, INBM |
401 | * |
402 | DO KZ=1,3 |
403 | VERTEX(KZ)=0 |
404 | ENDDO |
405 | IF(BEAMTY.EQ.1) THEN |
406 | CALL SXGAUS(RNDM,2) |
407 | VERTEX(1)= VERTEX(1)+FX+RNDM(1)*SIGX |
408 | VERTEX(2)= VERTEX(2)+FY+RNDM(2)*SIGY |
409 | VERTEX(3)= 0.05 |
410 | ELSEIF(BEAMTY.EQ.2) THEN |
411 | CALL GRNDM(RNDM,2) |
412 | VERTEX(1)= VERTEX(1)+FX+(2.*RNDM(1)-1.)*SIGX |
413 | VERTEX(2)= VERTEX(2)+FY+(2.*RNDM(2)-1.)*SIGY |
414 | VERTEX(3)= 0.05 |
415 | ENDIF |
416 | * |
417 | * Fill the vertex bank |
418 | * |
419 | CALL GSVERT(VERTEX,0,0,0.,0,NVTX) |
420 | IF(NVTX.EQ.0)WRITE(6,*)'Error defining vertex' |
421 | * |
422 | * Read HIJING event file |
423 | * |
424 | IF(IHIJ.GE.1)THEN |
425 | 10 READ(27,*,END=70) NMUL,EATT,JATT,NT,NP,N0,N01,N10,N11,B,NATT |
426 | IF(IHIJ.EQ.2)THEN |
427 | WRITE(6,*)'HIJING EVENT HEADER= ',NMUL,EATT,JATT,NT,NP, |
428 | + N0,N01,N10,N11,B,NATT |
429 | ENDIF |
430 | READ(27,*) (KATT(I,1),KATT(I,2),(PATT(I,J),J=1,4),I=1,NATT) |
431 | IVR=IVR+1 |
432 | IF(IVR.LT.IHIJF) GOTO 10 |
433 | * |
434 | BIMP=B |
435 | PARTI=FLOAT(NT)+FLOAT(NP) |
436 | RMULTI=FLOAT(NATT) |
437 | * |
438 | DO I=1,NATT |
439 | * -- Remove particles with negative z-direction |
440 | IF(PATT(I,3).LT.0.0) GOTO 30 |
441 | * -- Remove (if required) the spectators |
442 | IF(IHIJSP.EQ.1)THEN |
443 | IF(KATT(I,2).EQ.0..OR.KATT(I,2).EQ.10.) GOTO 30 |
444 | ENDIF |
445 | * -- Define participant flag : UBUF(1) |
446 | IF(KATT(I,2).EQ.0..OR.KATT(I,2).EQ.10.) THEN |
447 | UBUF(1)=0. |
448 | ELSE |
449 | UBUF(1)=1. |
450 | ENDIF |
451 | * -- Translate the particle ID from HIJING to GEANT |
452 | DO J=1,49 |
453 | IF(KATT(I,1).EQ.ITA(J)) THEN |
454 | IPART=J |
455 | GOTO 20 |
456 | ENDIF |
457 | ENDDO |
458 | 20 CONTINUE |
459 | * |
460 | DO JJ=1,3 |
461 | PLAB(JJ)=PATT(I,JJ) |
462 | ENDDO |
463 | * |
464 | * Rotations to account for beam divergence and crossing angle |
465 | IF(BMDIV.NE.0.) CALL ZDC_ROTP(PLAB,0) |
466 | IF(BMCRA.NE.0.) CALL ZDC_ROTP(PLAB,1) |
467 | * |
468 | * Fill the tracks bank |
469 | CALL GSKINE(PLAB,IPART,NVTX,UBUF,1,NT) |
470 | IF(IHIJ.EQ.2)THEN |
471 | WRITE(6,*)'IPART= ',IPART,'PLAB= ',PLAB |
472 | ENDIF |
473 | IF(NT.EQ.0)WRITE(6,*)'Error defining track' |
474 | * |
475 | 30 CONTINUE |
476 | ENDDO |
477 | * |
478 | * Read VENUS event file |
479 | * |
480 | ELSE IF(IVNU.GE.1)THEN |
481 | 40 READ(27,*,END=80)LEVT,NREVT,NMUL,B,KOLEVT,COLEVT, PMXEVT, |
482 | + EGYEVT,NP,NT,NATT |
483 | IF(LEVT.NE.1) THEN |
484 | WRITE(6,*) ' Error in reading VENUS data' |
485 | STOP |
486 | ENDIF |
487 | IF(IVNU.EQ.2)THEN |
488 | WRITE(6,*)'VENUS EVENT HEADER= ',LEVT,NREVT,NATT,B,KOLEVT, |
489 | + COLEVT,PMXEVT,EGYEVT,NP,NT |
490 | ENDIF |
491 | * |
492 | DO I=1,NATT |
493 | *-- read produced particles |
494 | READ(27,*,END=80) LPTL,NREVT,NRPTL,KATT(I,1), (PATT(I,J),J= |
495 | + 1,4),AMASSS |
496 | * |
497 | * lptl ......... record label (lptl=3) |
498 | * nrevt ........ event number |
499 | * nrptl ........ particle number |
500 | * katt(1,1)..... particle id |
501 | * patt ......... 4-momentum (px,py,pz,en) in lab |
502 | * amasss ........ mass of particle |
503 | * |
504 | IF(LPTL.NE.3) THEN |
505 | WRITE(6,*) ' ERROR IN READING VENUS DATA' |
506 | STOP |
507 | ENDIF |
508 | ENDDO |
509 | * |
510 | IVR=IVR+1 |
511 | IF(IVR.LT.IVNUF) GOTO 40 |
512 | * |
513 | BIMP=B |
514 | PARTI=(FLOAT(NT)+FLOAT(NP))/2. |
515 | RMULTI=FLOAT(NATT) |
516 | * |
517 | DO I=1,NATT |
518 | * -- remove particles in negative direction |
519 | IF(PATT(I,3).LT.0.0) GOTO 60 |
520 | * -- if required remove the spectators |
521 | IF(IVNSP.EQ.1)THEN |
522 | IF((KATT(I,1).EQ.1120.OR.KATT(I,1).EQ.1220).AND. (PATT(I |
523 | + ,4).LT.(EPERN*1.005).AND. PATT(I,4).GT.(EPERN*0.995))) |
524 | + GOTO 60 |
525 | ENDIF |
526 | * -- flag the spectators (UBUF(1)=0) |
527 | IF((KATT(I,1).EQ.1120.OR.KATT(I,1).EQ.1220).AND. (PATT(I,4) |
528 | + .GT.(EPERN*0.995).AND. PATT(I,4).LT.(EPERN*1.005))) THEN |
529 | UBUF(1)=0. |
530 | ELSE |
531 | UBUF(1)=1. |
532 | ENDIF |
533 | * -- Translate the particle ID from VENUS to GEANT |
534 | DO J=1,49 |
535 | IF(KATT(I,1).EQ.ITB(J))THEN |
536 | IPART=J |
537 | GOTO 50 |
538 | ENDIF |
539 | ENDDO |
540 | 50 CONTINUE |
541 | DO JJ=1,3 |
542 | PLAB(JJ)=PATT(I,JJ) |
543 | ENDDO |
544 | * |
545 | * Rotations to account for beam divergence and crossing angle |
546 | IF(BMCRA.NE.0.)CALL ZDC_ROTP(PLAB,1) |
547 | IF(BMDIV.NE.0.)CALL ZDC_ROTP(PLAB,0) |
548 | * |
549 | * Fill the tracks bank |
550 | CALL GSKINE(PLAB,IPART,NVTX,UBUF,1,NT) |
551 | IF(IVNU.EQ.2)THEN |
552 | WRITE(6,*)'IPART= ',IPART,'PLAB= ',PLAB |
553 | ENDIF |
554 | IF(NT.EQ.0)WRITE(6,*)'Error defining track' |
555 | 60 CONTINUE |
556 | ENDDO |
557 | * |
558 | ELSE |
559 | * |
560 | * Generation of test particles |
561 | * |
562 | * Definition of the incident particle |
563 | * |
564 | IPART = NINT(PKINE(1)) |
565 | UBUF(1)=1 |
566 | CALL GFPART(IPART,PARTNAME,ITRTYP,AMASS,CHARGE,TLIFE,UB,NB) |
567 | AAP=1.0 |
568 | * |
569 | * --- In case of a ion , find out energy and momentum |
570 | * --- NAAP=A=Z+N - atomic number of the initial particle |
571 | IF(IPART.GE.45) THEN |
572 | NAAP = INT(0.5 + AMASS/0.93149432) |
573 | AAP = NAAP |
574 | UBUF(1)=0. |
575 | ENDIF |
576 | IF((IPART.GE.45).OR.(IPART.EQ.13).OR.(IPART.EQ.14)) THEN |
577 | UBUF(1)=0. |
578 | ENDIF |
579 | * |
580 | IF(IPART.EQ.14)THEN |
581 | AMA=0.938 |
582 | ELSE IF(IPART.EQ.13)THEN |
583 | AMA=0.939 |
584 | ELSE |
585 | AMA=AMASS |
586 | ENDIF |
587 | * |
588 | IF(PKINE(6).EQ.0.)THEN |
589 | PLAB(1)=PKINE(2)*PKINE(3) |
590 | PLAB(2)=PKINE(2)*PKINE(4) |
591 | PLAB(3)=PKINE(2)*PKINE(5) |
592 | ELSE |
593 | SCANG=2*ATAN(EXP(-(PKINE(6)))) |
594 | PLAB(1)=-PKINE(2)*SIN(SCANG) |
595 | PLAB(2)=0. |
596 | PLAB(3)=PKINE(2)*COS(SCANG) |
597 | ENDIF |
598 | * |
599 | * Rotations to account for beam divergence and crossing angle |
600 | IF(BMCRA.NE.0.)CALL ZDC_ROTP(PLAB,1) |
601 | IF(BMDIV.NE.0.)CALL ZDC_ROTP(PLAB,0) |
602 | * |
603 | * -- in case of a nucleon, if required, apply the Fermi momentum |
604 | IF(UBUF(1).EQ.0.AND.IFERMI.EQ.1) THEN |
605 | IF(IPART.EQ.13.OR.IPART.EQ.14)THEN |
606 | CALL ZDC_EFERMI(IPART) |
607 | ELSE IF(IPART.GE.45)THEN |
608 | CALL ZDC_EFERMI(13) |
609 | GOLDHABER=SQRT(AAP*(207-AAP)/(207-1.)) |
610 | DO JJ=1,3 |
611 | DDP(JJ)=DDP(JJ)*GOLDHABER |
612 | ENDDO |
613 | ENDIF |
614 | DO II=1,3 |
615 | BALP(II)=-PLAB(II) |
616 | ENDDO |
617 | BALP(4)=SQRT(BALP(1)**2+BALP(2)**2+BALP(3)**2+AMA**2) |
618 | DO II=1,3 |
619 | DDDP(II)=DDP(II) |
620 | ENDDO |
621 | DDDP(4)=SQRT(DDDP(1)**2+DDDP(2)**2+DDDP(3)**2+AMA**2) |
622 | CALL LORENF(AMA,BALP,DDDP,PLABS) |
623 | DO II=1,3 |
624 | PLAB(II)=PLABS(II) |
625 | ENDDO |
626 | c WRITE(6,*)'BALP ',BALP,'DDDP ',DDDP,'PLABS ',PLABS |
627 | ENDIF |
628 | * |
629 | BIMP=999. |
630 | * |
631 | * Fill the tracks bank |
632 | CALL GSKINE(PLAB,IPART,NVTX,UBUF,1,NT) |
633 | IF(NT.EQ.0)WRITE(6,*)'Error defining track' |
634 | |
635 | ENDIF |
636 | * |
637 | * Some debug printings |
638 | IF(IHIJ.EQ.2.OR.IVNU.EQ.2)THEN |
639 | CALL GPVERT(NVTX) |
640 | DO KK=1,NT |
641 | CALL GPKINE(KK) |
642 | ENDDO |
643 | ENDIF |
644 | * |
645 | C |
646 | RETURN |
647 | 70 WRITE(6,*) 'END OF HIJING FILE' |
648 | IEORUN=1 |
649 | GOTO 999 |
650 | 80 WRITE(6,*) 'END OF VENUS FILE' |
651 | IEORUN=1 |
652 | * |
653 | 999 END |
654 | |
655 | *CMZ : 19/11/98 09.54.52 by Federico Carminati |
656 | *-- Author : |
657 | SUBROUTINE ZDC_ROTP(PLAB,ICROSS) |
658 | #undef CERNLIB_GEANT321_GCKINE_INC |
659 | #include "geant321/gckine.inc" |
660 | #undef CERNLIB_GEANT321_GCFLAG_INC |
661 | #include "geant321/gcflag.inc" |
662 | C |
663 | #undef CERNLIB_GEANT321_GCONSP_INC |
664 | #include "geant321/gconsp.inc" |
665 | #include "zdc_common.inc" |
666 | C |
667 | DOUBLE PRECISION TETPART,FIPART,TETDIV,FIDIV,TETSUM,FISUM |
668 | DOUBLE PRECISION DPLAB(3) |
669 | DIMENSION PLAB(3) |
670 | DIMENSION RVEC(1) |
671 | DIMENSION RNDM(1) |
672 | C |
673 | C TAKE NON-ZERO DIVERGENCY |
674 | c write(6,*)'PLAB (NODIV) = ',plab |
675 | PMOD=0.0 |
676 | DO IJ=1,3 |
677 | DPLAB(IJ)=PLAB(IJ) |
678 | PMOD=PMOD+PLAB(IJ)*PLAB(IJ) |
679 | ENDDO |
680 | PMOD=SQRT(PMOD) |
681 | IF(ICROSS.EQ.0)THEN |
682 | C--- TETDIV is given in radiants and it is read by data cards |
683 | CALL RNORML(RVEC,1) |
684 | c WRITE(6,*)RVEC |
685 | TETDIV=PKINE(10)*ABS(RVEC(1)) |
686 | C--- PHIDIV is randomly chosen between 0 and 360 |
687 | CALL GRNDM(RNDM,1) |
688 | FIDIV=RNDM(1)*TWOPI |
689 | ELSE IF(ICROSS.EQ.1)THEN |
690 | IF(ISWIT(7).EQ.0)THEN |
691 | TETDIV=0. |
692 | FIDIV=0. |
693 | ELSE IF(ISWIT(7).EQ.1)THEN |
694 | TETDIV=PKINE(4) |
695 | FIDIV=0. |
696 | ELSE IF(ISWIT(7).EQ.2)THEN |
697 | TETDIV=PKINE(4) |
698 | FIDIV=PIBY2 |
699 | ENDIF |
700 | ENDIF |
701 | TETPART=DATAN(DSQRT(DPLAB(1)**2+DPLAB(2)**2)/DPLAB(3)) |
702 | IF(DPLAB(2).NE.0..OR.DPLAB(1).NE.0.) THEN |
703 | FIPART=DATAN2(DPLAB(2),DPLAB(1)) |
704 | ELSE |
705 | FIPART=0. |
706 | ENDIF |
707 | IF (FIPART.LT.0)FIPART=FIPART+TWOPI |
708 | C write(6,*)' TETPART= ',tetpart,' FIPART= ',fipart |
709 | TETDIV=TETDIV*RADDEG |
710 | FIDIV=FIDIV*RADDEG |
711 | TETPART=TETPART*RADDEG |
712 | FIPART=FIPART*RADDEG |
713 | c CALL ZDC_ADDANG(TETDIV,FIDIV,TETPART,FIPART,TETSUM,FISUM) |
714 | CALL ZDC_ADDANG(TETPART,FIPART,TETDIV,FIDIV,TETSUM,FISUM) |
715 | C write(6,*)' TETSUM= ',tetsum,'FISUM= ',fisum |
716 | TETSUM=TETSUM*DEGRAD |
717 | FISUM=FISUM*DEGRAD |
718 | PLAB(1)=PMOD*DSIN(TETSUM)*DCOS(FISUM) |
719 | PLAB(2)=PMOD*DSIN(TETSUM)*DSIN(FISUM) |
720 | PLAB(3)=PMOD*DCOS(TETSUM) |
721 | C write(6,*)'PLAB (DIV) = ',plab |
722 | RETURN |
723 | END |
724 | |
725 | *CMZ : 19/11/98 09.54.52 by Federico Carminati |
726 | *CMZ : 2.03/01 25/08/98 23.40.04 by Federico Carminati |
727 | *-- Author : E.SCOMPARIN, 13/05/1996. |
728 | SUBROUTINE ZDC_STEP |
729 | C |
730 | C *** RECORDING OF ZDC HITS AFTER EACH TRACKING STEP *** |
731 | C |
732 | C CALLED BY : SXSTEP |
733 | C ORIGIN : E. SCOMPARIN |
734 | C |
735 | #undef CERNLIB_GEANT321_GCKINE_INC |
736 | #undef CERNLIB_GEANT321_GCONSP_INC |
737 | #undef CERNLIB_GEANT321_GCVOLU_INC |
738 | #undef CERNLIB_GEANT321_GCKMAX_INC |
739 | #undef CERNLIB_GEANT321_GCKING_INC |
740 | #undef CERNLIB_GEANT321_GCNUM_INC |
741 | #undef CERNLIB_GEANT321_GCTRAK_INC |
742 | #undef CERNLIB_GEANT321_GCFLAG_INC |
743 | #undef CERNLIB_GEANT321_GCTMED_INC |
744 | #undef CERNLIB_GEANT321_GCMATE_INC |
745 | #undef CERNLIB_GEANT321_GCSETS_INC |
746 | #include "geant321/gckine.inc" |
747 | #include "geant321/gconsp.inc" |
748 | #include "geant321/gcvolu.inc" |
749 | #include "geant321/gcking.inc" |
750 | #include "geant321/gcnum.inc" |
751 | #include "geant321/gctrak.inc" |
752 | #include "geant321/gcflag.inc" |
753 | #include "geant321/gctmed.inc" |
754 | #include "geant321/gcmate.inc" |
755 | #include "geant321/gcsets.inc" |
756 | #include "zdc_common.inc" |
757 | * |
758 | DIMENSION LNAM(15),LNUM(15) |
759 | DIMENSION UB(1) |
760 | * |
761 | DIMENSION VERT0(3),PVERT0(4) |
762 | CHARACTER*20 CHNPART0 |
763 | * |
764 | DIMENSION XLOCM(3),XANGM(3),XLOCL(3),XANGL(3),FNULA(3),ZENTR(3) |
765 | DATA ALFMAX /110.0/ |
766 | * |
767 | * New track |
768 | * |
769 | IF(ITRA.NE.ITOLD)THEN |
770 | IF(INWVOL.EQ.1.AND.IDTYPE.GT.0)THEN |
771 | * |
772 | * Impact point on the ZDC (either ZPRO or ZNEU) |
773 | * |
774 | NLEV0=NLEVEL |
775 | CALL GFPATH(ISET,IDET,NUMBV,NLEV,LNAM,LNUM) |
776 | CALL UCTOH('ZNEU',IZNEU,4,4) |
777 | CALL UCTOH('ZPRO',IZPRO,4,4) |
778 | DO I=NLEV,1,-1 |
779 | IF(LNAM(I).EQ.IZNEU.OR.LNAM(I).EQ.IZPRO)THEN |
780 | NLEVEL=0 |
781 | CALL GLVOLU(I,LNAM,LNUM,IER) |
782 | IF(IDTYPE.EQ.1)THEN |
783 | CALL GMTOD(VECT,HITSBUPN(3),1) |
784 | HITSBUPP(3)=99. |
785 | HITSBUPP(4)=99. |
786 | ELSE IF(IDTYPE.EQ.2)THEN |
787 | CALL GMTOD(VECT,HITSBUPP(3),1) |
788 | HITSBUPN(3)=99. |
789 | HITSBUPN(4)=99. |
790 | ENDIF |
791 | NLEVEL=NLEV0 |
792 | GOTO 10 |
793 | ENDIF |
794 | ENDDO |
795 | 10 CONTINUE |
796 | * |
797 | * Particle type and energy for the primary track |
798 | * |
799 | CALL GFKINE(ITRA,VERT0,PVERT0,IPART0,NVERT0,UB,NB) |
800 | CALL GFPART(IPART0,CHNPART0,ITRTYP0,AMASS0,CHARGE0,TLIFE0, |
801 | + UB,NB) |
802 | PVERTM=VMOD(PVERT0,3) |
803 | GEKIN0=SQRT(PVERTM**2+AMASS0**2) |
804 | IF(IPART.NE.IPART0.OR.GEKIN.LT.GEKIN0*0.95)THEN |
805 | HITSBUPP(1)=FLOAT(IPART0) |
806 | HITSBUPP(2)=GEKIN0 |
807 | HITSBUPP(5)=1. |
808 | HITSBUPN(1)=FLOAT(IPART0) |
809 | HITSBUPN(2)=GEKIN0 |
810 | HITSBUPN(5)=1. |
811 | ELSE |
812 | HITSBUPP(1)=FLOAT(IPART) |
813 | HITSBUPP(2)=GEKIN |
814 | HITSBUPP(5)=0. |
815 | HITSBUPN(1)=FLOAT(IPART) |
816 | HITSBUPN(2)=GEKIN |
817 | HITSBUPN(5)=0. |
818 | ENDIF |
819 | HITSBUPP(6)=0. |
820 | HITSBUPP(7)=0. |
821 | HITSBUPN(6)=0. |
822 | HITSBUPN(7)=0. |
823 | IF(IDTYPE.EQ.1)THEN |
824 | CALL GSAHIT(ISET,IDET,ITRA,NUMBV,HITSBUPN,IHIT) |
825 | ELSE IF(IDTYPE.EQ.2)THEN |
826 | CALL GSAHIT(ISET,IDET,ITRA,NUMBV,HITSBUPP,IHIT) |
827 | ENDIF |
828 | ITOLD=ITRA |
829 | ENDIF |
830 | ENDIF |
831 | * |
832 | * Energy deposition in the passive material |
833 | * |
834 | IF(INWVOL.EQ.0.AND.IDTYPE.GT.0)THEN |
835 | HITS(6)=0. |
836 | IF(ISTOP.EQ.2)THEN |
837 | HITS(7)=GEKIN |
838 | ELSE |
839 | HITS(7)=DESTEP |
840 | ENDIF |
841 | IF(IDTYPE.EQ.1)THEN |
842 | DO KCOPY=1,5 |
843 | HITSBUPN(KCOPY)=HITS(KCOPY) |
844 | ENDDO |
845 | ELSE IF(IDTYPE.EQ.2)THEN |
846 | DO KCOPY=1,5 |
847 | HITSBUPP(KCOPY)=HITS(KCOPY) |
848 | ENDDO |
849 | ENDIF |
850 | CALL GSCHIT(ISET,IDET,ITRA,NUMBV,HITS,NHSUM,IHIT) |
851 | ENDIF |
852 | * |
853 | * Select charged relativistic particles |
854 | * |
855 | IF(ABS(CHARGE).EQ.0)RETURN |
856 | IF(GEKIN.GT.0.00001) THEN |
857 | BETA = VECT(7)/GETOT |
858 | ELSE |
859 | RETURN |
860 | ENDIF |
861 | IF(BETA.LT.0.67) RETURN |
862 | |
863 | * |
864 | * Particle inside fiber |
865 | * |
866 | IF(NUMED.EQ.803.AND.IDTYPE.GT.0) THEN |
867 | * |
868 | XLOCM(1)= VECT(1) |
869 | XLOCM(2)= VECT(2) |
870 | XLOCM(3)= VECT(3) |
871 | XANGM(1)= VECT(4) |
872 | XANGM(2)= VECT(5) |
873 | XANGM(3)= VECT(6) |
874 | * |
875 | * --- Calculating the position of the centre of the fiber in MARS |
876 | DO KZERO=1,3 |
877 | FNULA(KZERO)=0 |
878 | ZENTR(KZERO)=0 |
879 | ENDDO |
880 | CALL GDTOM(FNULA,ZENTR,1) |
881 | * |
882 | * --- Angle between the particle trajectory and the fiber axis |
883 | * |
884 | CALL GMTOD(XANGM,XANGL,2) |
885 | IF(XANGL(3).GT.1.0) XANGL(3)=SIGN(1.0,XANGL(3)) |
886 | ALFAR=ACOS(XANGL(3)) |
887 | ALFAD=ALFAR*RADDEG |
888 | * |
889 | IF (ALFAD.GT.ALFMAX) RETURN |
890 | IALFA = INT(1.+ALFAD/2.) |
891 | * |
892 | * --- Distance between the particle trajectory and the fiber axis |
893 | * |
894 | CALL GMTOD(XLOCM,XLOCL,1) |
895 | IF(ABS(XANGL(1)).GT.0.0001) THEN |
896 | DCOEFF = XANGL(2)/XANGL(1) |
897 | BE=ABS(XLOCL(2)-DCOEFF*XLOCL(1))/SQRT(DCOEFF*DCOEFF+1.) |
898 | ELSE |
899 | BE = ABS(XLOCL(1)) |
900 | ENDIF |
901 | IF(IDTYPE.EQ.1)THEN |
902 | RADIUS=FIZN(2) |
903 | ELSEIF(IDTYPE.EQ.2)THEN |
904 | RADIUS=FIZP(2) |
905 | ENDIF |
906 | IF(BE.GT.RADIUS) THEN |
907 | PRINT*,'DISTANCE GREATER THAN FIBER RADIUS' |
908 | RETURN |
909 | ENDIF |
910 | IBE = 1+INT(1000.*BE) |
911 | IF(IBE.GT.NBE)IBE=NBE |
912 | * |
913 | * --- Finding out proper index for beta |
914 | IF(BETA.GE.0.67 .AND. BETA.LE.0.75) IBETA=1 |
915 | IF(BETA.GT.0.75 .AND. BETA.LE.0.85) IBETA=2 |
916 | IF(BETA.GT.0.85 .AND. BETA.LE.0.95) IBETA=3 |
917 | IF(BETA.GT.0.95 .AND. BETA.LE.1.00) IBETA=4 |
918 | * |
919 | * --- Looking into the tables |
920 | OUT = CHARGE*CHARGE*TABLE(IBETA,IALFA,IBE) |
921 | CALL RNPSSN(OUT,N,IERROR) |
922 | HITS(6)=FLOAT(N) |
923 | HITS(7)=0. |
924 | IF(IDTYPE.EQ.1)THEN |
925 | CALL UCOPY(HITSBUPN,HITS,4) |
926 | ELSE IF(IDTYPE.EQ.2)THEN |
927 | CALL UCOPY(HITSBUPP,HITS,4) |
928 | ENDIF |
929 | CALL GSCHIT(ISET,IDET,ITRA,NUMBV,HITS,NHSUM,IHIT) |
930 | ENDIF |
931 | * |
932 | 999 END |
933 | |
934 | *CMZ : 18/11/98 16.55.53 by Federico Carminati |
935 | *CMZ : 2.03/01 06/08/98 10.13.59 by Federico Carminati |
936 | *-- Author : E. SCOMPARIN |
937 | SUBROUTINE ZDC_TRKI(IFLAG) |
938 | C |
939 | C *** DEFINITION OF THE GEOMETRY OF THE PHOS *** |
940 | C *** NVE 24-SEP-1990 CERN GENEVA *** |
941 | C |
942 | C CALLED BY : SXTRKI |
943 | C ORIGIN : E. SCOMPARIN |
944 | C |
945 | #undef CERNLIB_GEANT321_GCKINE_INC |
946 | #undef CERNLIB_GEANT321_GCKMAX_INC |
947 | #undef CERNLIB_GEANT321_GCKING_INC |
948 | #undef CERNLIB_GEANT321_GCTRAK_INC |
949 | #undef CERNLIB_GEANT321_GCTMED_INC |
950 | #undef CERNLIB_GEANT321_GCFLAG_INC |
951 | #undef CERNLIB_GEANT321_GCMATE_INC |
952 | #include "geant321/gckine.inc" |
953 | #include "geant321/gcking.inc" |
954 | #include "geant321/gctrak.inc" |
955 | #include "geant321/gctmed.inc" |
956 | #include "geant321/gcflag.inc" |
957 | #include "geant321/gcmate.inc" |
958 | #include "zdc_common.inc" |
959 | |
960 | REAL UBUF(1) |
961 | DATA ITOLDG/0/ |
962 | |
963 | * Fetch new track (VENUS or HIJING) |
964 | IF(ITRA.NE.ITOLDG)THEN |
965 | CALL GFKINE(ITRA,VERT,PVERT,IPART,NVERT,UBUF,NUBUF) |
966 | IFLPART=INT(UBUF(1)) |
967 | IF(IHIJ.EQ.2.OR.IVNU.EQ.2)WRITE(6,*)'ITRA,IFLPART ',ITRA, |
968 | + IFLPART |
969 | ITOLDG=ITRA |
970 | ENDIF |
971 | * |
972 | END |
973 | |
974 | *CMZ : 15/02/99 14.40.30 by Federico Carminati |
975 | *-- Author : Rene Brun |
976 | SUBROUTINE ZDC_DEFS |
977 | C |
978 | C *** DEFINITION OF THE ZDC parameters |
979 | C CALLED BY : SXDEFS |
980 | C |
981 | end |
982 | |
983 | subroutine zdc_setbeam(beam, fxoff, fyoff, sx, sy, |
984 | + div, angle, cross) |
985 | * beam : 1 = gaussian beam |
986 | * : 2 = uniform beam |
987 | * fxoff : x-coordinate of beam offset |
988 | * fyoff : y-coordinate of beam offset |
989 | * sx : sigma-x of the beam (gaussian or uniform) |
990 | * sy : sigma-y of the beam (gaussian or uniform) |
991 | * div : divergency of the beam (32*10**-6 rad for LHC) |
992 | * angle : beam crossing angle (100*10**-6 rad for LHC) |
993 | * cross : 1 = horizontal beam crossing |
994 | * : 2 = vertical beam crossing |
995 | #include "zdc_common.inc" |
996 | integer beam, cross |
997 | real fxoff,fyoff,sx,sy,div,angle |
998 | beamty = type |
999 | fx = fxoff |
1000 | fy = fyoff |
1001 | sigx = sx |
1002 | sigy = sy |
1003 | bmdiv = div |
1004 | bmcra = angle |
1005 | iflcr = cross |
1006 | end |
1007 | |
1008 | subroutine zdc_sethijing(hij, hijf, hijsp, file) |
1009 | * IHIJ : 1 = read HIJING event file |
1010 | * : 2 = " " " " + debug |
1011 | * IHIJF : event number of the first event to be read from file |
1012 | * IHIJSP: 0 = read all particles |
1013 | * : 1 = remove spectator nucleons |
1014 | #include "zdc_common.inc" |
1015 | integer hij,hijf,hijsp |
1016 | character *(*) file |
1017 | ihij = hij |
1018 | ihijf = hijf |
1019 | ihijsp = hijsp |
1020 | fileh = file |
1021 | end |
1022 | |
1023 | subroutine zdc_setvenus(hiv, hivf, hivsp, file) |
1024 | * IHIV : 1 = read VENUS event file |
1025 | * : 2 = " " " " + debug |
1026 | * IHIVF : event number of the first event to be read from file |
1027 | * IHIVSP: 0 = read all particles |
1028 | * : 1 = remove spectator nucleons |
1029 | #include "zdc_common.inc" |
1030 | integer hiv,hivf,hivsp |
1031 | character *(*) file |
1032 | ihiv = hiv |
1033 | ihivf = hivf |
1034 | ihivsp = hivsp |
1035 | filev = file |
1036 | end |
1037 | |
1038 | subroutine zdc_setkine(code, pmom, cx, cy, cz, type, fermi) |
1039 | * code : GEANT code of the test particle |
1040 | * pmom : absolute value of particle momentum |
1041 | * cx,cy,cz : director cosines of the track (if type) |
1042 | * type : 0 = take director cosines from cx,cy,cz) |
1043 | * : <>0 = pseudorapidity of the test particle |
1044 | * fermi : 0 = no Fermi motion for the spectator nucleons |
1045 | * : 1 = Fermi motion for the spectator nucleons |
1046 | #undef CERNLIB_GEANT321_GCKINE_INC |
1047 | #include "geant321/gckine.inc" |
1048 | #include "zdc_common.inc" |
1049 | integer code,type,fermi |
1050 | real pmom,cx,cy,cz |
1051 | pkine(1) = code |
1052 | pkine(2) = pmom |
1053 | pkine(3) = cx |
1054 | pkine(4) = cy |
1055 | pkine(5) = cz |
1056 | pkine(6) = type |
1057 | ifermi = fermi |
1058 | end |
1059 | *CMZ : 19/11/98 09.31.36 by Federico Carminati |
1060 | *-- Author : Federico Carminati 19/11/98 |
1061 | SUBROUTINE SXGAUS(RVEC,N) |
1062 | #undef CERNLIB_GEANT321_GCONSP_INC |
1063 | #include "geant321/gconsp.inc" |
1064 | *KEND. |
1065 | * |
1066 | * Vector of random numbers distributed like a gaussian |
1067 | DIMENSION RNDM(2), RVEC(N) |
1068 | DO K=1,(N-1)/2+1 |
1069 | CALL GRNDM(RNDM,2) |
1070 | PHI = TWOPI*RNDM(1) |
1071 | RHO = SQRT(-2.*LOG(RNDM(2))) |
1072 | RVEC(2*K-1) = RHO*SIN(PHI) |
1073 | IF(2*K.LE.N) RVEC(2*K ) = RHO*COS(PHI) |
1074 | ENDDO |
1075 | * |
1076 | END |