]> git.uio.no Git - u/mrichter/AliRoot.git/blame - HERWIG/src/hwhppe.f
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Coding rule violations corrected.
[u/mrichter/AliRoot.git] / HERWIG / src / hwhppe.f
CommitLineData
3820ca8e 1
2CDECK ID>, HWHPPE.
3
4*CMZ :- -12/01/93 10.12.43 by Bryan Webber
5
6*-- Author : Ian Knowles
7
8C-----------------------------------------------------------------------
9
10 SUBROUTINE HWHPPE
11
12C-----------------------------------------------------------------------
13
14C point-like photon/QCD heavy flavour single excitation, using exact
15
16C massive lightcone kinematics, mean EVWGT = sigma in nb.
17
18C-----------------------------------------------------------------------
19
20 INCLUDE 'HERWIG61.INC'
21
22 DOUBLE PRECISION HWR,HWRUNI,HWUALF,EPS,PP1,PP2,QM2,FACTR,
23
24 & PT,PJ,PT2,PTM,EXY,T,CC,EXY2,S,U,C,SIGE,HCS,RCS
25
26 INTEGER IQ1,IQ2,ID1,ID2,IHAD1,IHAD2
27
28 EXTERNAL HWR,HWRUNI,HWUALF
29
30 SAVE PP1,PP2,IQ1,IQ2,QM2,FACTR,SIGE,HCS
31
32 PARAMETER (EPS=1.E-9)
33
34 IHAD1=1
35
36 IF (JDAHEP(1,IHAD1).NE.0) IHAD1=JDAHEP(1,IHAD1)
37
38 IHAD2=2
39
40 IF (JDAHEP(1,IHAD2).NE.0) IHAD2=JDAHEP(1,IHAD2)
41
42 IF (FSTWGT.OR.IHAD1.NE.1.OR.IHAD2.NE.2) THEN
43
44 PP1=PHEP(4,IHAD1)+ABS(PHEP(3,IHAD1))
45
46 PP2=PHEP(4,IHAD2)+ABS(PHEP(3,IHAD2))
47
48 XX(1)=1.
49
50 IQ1=MOD(IPROC,100)
51
52 IQ2=IQ1+6
53
54 QM2=RMASS(IQ1)**2
55
56 FACTR=GEV2NB*(YJMAX-YJMIN)*4.*PIFAC*CFFAC*PP1
57
58 & *ALPHEM*QFCH(IQ1)**2
59
60 ENDIF
61
62 IF (GENEV) THEN
63
64 RCS=HCS*HWR()
65
66 ELSE
67
68 EVWGT=0.
69
70 CALL HWRPOW(PT,PJ)
71
72 PT2=PT**2
73
74 PTM=SQRT(PT2+QM2)
75
76 EXY=EXP(HWRUNI(1,YJMIN,YJMAX))
77
78 T=-PP1*PT/EXY
79
80 CC=T**2-4.*QM2*(PT2+T)
81
82 IF (CC.LT.ZERO) RETURN
83
84 EXY2=(2.*PT2+T-SQRT(CC))*PP1/(2.*T*PTM)
85
86 IF (EXY2.LT.EXP(YJMIN).OR.EXY2.GT.EXP(YJMAX)) RETURN
87
88 XX(2)=(PT/EXY+PTM/EXY2)/PP2
89
90 IF (XX(2).GT.ONE) RETURN
91
92C define: S=Shat-M**2, T=That ,U=Uhat-M**2 (2p.Q, -2p.g, -2p.Q')
93
94 S=XX(2)*PP1*PP2
95
96 U=-S-T
97
98 COSTH=(1.+QM2/S)*(T-U)/S-QM2/S
99
100C Set hard process scale (Approx ET-jet)
101
102 EMSCA=SQRT(2.*S*T*U/(S*S+T*T+U*U))
103
104 C=QM2*T/(U*S)
105
106 SIGE=-FACTR*PT*PJ*HWUALF(1,EMSCA)*(S/U+U/S+4.*C*(1.-C))
107
108 & /(S**2*EXY2*PTM*(1-QM2/(XX(2)*PP2*EXY2)**2))
109
110 CALL HWSFUN(XX(2),EMSCA,IDHW(IHAD2),NSTRU,DISF(1,2),2)
111
112 ENDIF
113
114 HCS=0.
115
116 ID1=59
117
118C photon+Q ---> g+Q
119
120 ID2=IQ1
121
122 IF (DISF(ID2,2).LT.EPS) GOTO 10
123
124 HCS=HCS+SIGE*DISF(ID2,2)
125
126 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(13,ID2,1423,51,*99)
127
128C photon+Qbar ---> g+Qbar
129
130 10 ID2=IQ2
131
132 IF (DISF(ID2,2).LT.EPS) GOTO 20
133
134 HCS=HCS+SIGE*DISF(ID2,2)
135
136 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(13,ID2,1342,52,*99)
137
138 20 EVWGT=HCS
139
140 RETURN
141
142C Generate event
143
144 99 IDN(1)=ID1
145
146 IDN(2)=ID2
147
148 IDCMF=15
149
150 CALL HWETWO
151
152 999 END
153
154CDECK ID>, HWHPPH.
155
156*CMZ :- -12/01/93 10.12.43 by Bryan Webber
157
158*-- Author : Ian Knowles
159
160C-----------------------------------------------------------------------
161
162 SUBROUTINE HWHPPH
163
164C-----------------------------------------------------------------------
165
166C Point-like photon/gluon heavy flavour pair production, with
167
168C exact lightcone massive kinematics, mean EVWGT = sigma in nb.
169
170C-----------------------------------------------------------------------
171
172 INCLUDE 'HERWIG61.INC'
173
174 DOUBLE PRECISION HWRUNI,HWUALF,EPS,PP1,PP2,QM2,FACTR,ET,EJ,ET2,
175
176 & EXY,EXY2,S,T,U,C
177
178 INTEGER IQ1,IHAD1,IHAD2
179
180 EXTERNAL HWRUNI,HWUALF
181
182 SAVE PP1,PP2,IQ1,QM2,FACTR
183
184 PARAMETER (EPS=1.E-9)
185
186 IHAD1=1
187
188 IF (JDAHEP(1,IHAD1).NE.0) IHAD1=JDAHEP(1,IHAD1)
189
190 IHAD2=2
191
192 IF (JDAHEP(1,IHAD2).NE.0) IHAD2=JDAHEP(1,IHAD2)
193
194 IF (FSTWGT.OR.IHAD1.NE.1.OR.IHAD2.NE.2) THEN
195
196 PP1=PHEP(4,IHAD1)+ABS(PHEP(3,IHAD1))
197
198 PP2=PHEP(4,IHAD2)+ABS(PHEP(3,IHAD2))
199
200 XX(1)=1.
201
202 IQ1=MOD(IPROC,100)
203
204 QM2=RMASS(IQ1)**2
205
206 IHPRO=53
207
208 FACTR=-GEV2NB*(YJMAX-YJMIN)*.5*PIFAC*ALPHEM*QFCH(IQ1)**2
209
210 ENDIF
211
212 IF (GENEV) THEN
213
214C Generate event
215
216 IDN(1)=59
217
218 IDN(2)=13
219
220 IDN(3)=IQ1
221
222 IDN(4)=IQ1+6
223
224 ICO(1)=1
225
226 ICO(2)=4
227
228 ICO(3)=2
229
230 ICO(4)=3
231
232 IDCMF=15
233
234 CALL HWETWO
235
236 ELSE
237
238C Select kinematics
239
240 EVWGT=0.
241
242 CALL HWRPOW(ET,EJ)
243
244 ET2=ET**2
245
246 EXY=EXP(HWRUNI(1,YJMIN,YJMAX))
247
248 EXY2=2.*PP1/ET-EXY
249
250 IF (EXY2.LT.EXP(YJMIN).OR.EXY2.GT.EXP(YJMAX)) RETURN
251
252 XX(2)=.5*ET*(1./EXY+1./EXY2)/PP2
253
254 IF (XX(2).LT.ZERO.OR.XX(2).GT.ONE) RETURN
255
256 S=XX(2)*PP1*PP2
257
258 IF (S.LT.ET2) RETURN
259
260C define: S=Shat, T=That-M**2, U=Uhat-M**2 (2p.g, -2p.Q, -2p.QBar)
261
262 T=-.5*PP1*ET/EXY
263
264 U=-S-T
265
266 COSTH=(T-U)/(S*SQRT(1.-4.*QM2/S))
267
268 EMSCA=SQRT(2.*S*T*U/(S*S+T*T+U*U))
269
270 CALL HWSFUN(XX(2),EMSCA,IDHW(IHAD2),NSTRU,DISF(1,2),2)
271
272C photon+g ---> Q+Qbar
273
274 IF (DISF(13,2).LT.EPS) THEN
275
276 EVWGT=0.
277
278 ELSE
279
280 C=QM2*S/(U*T)
281
282 EVWGT=FACTR*EJ*ET*HWUALF(1,EMSCA)
283
284 & *DISF(13,2)*(T/U+U/T+4.*C*(1.-C))/(S*T)
285
286 ENDIF
287
288 ENDIF
289
290 999 END
291
292CDECK ID>, HWHPPM.
293
294*CMZ :- -09/12/93 15.50.26 by Mike Seymour
295
296*-- Author : Ian Knowles & Mike Seymour
297
298C-----------------------------------------------------------------------
299
300 SUBROUTINE HWHPPM
301
302C-----------------------------------------------------------------------
303
304C Point-like photon/QCD direct meson production
305
306C See M. Benayoun, et al., Nucl. Phys. B282 (1987) 653 for details.
307
308C mean EVWGT = sigma in nb
309
310C-----------------------------------------------------------------------
311
312 INCLUDE 'HERWIG61.INC'
313
314 DOUBLE PRECISION HWR,HWRUNI,HWUALF,EPS,PP1,PP2,ET,EJ,EXY,EXY2,
315
316 & FACT,FACTR,S,T,U,REDS,DELT(3,3),C1STU,C3STU,HCS,RCS,CMIX,SMIX,
317
318 & C1WVFN,FPI,FETA8,FETA1,FRHO,FPHI8,FPHI1,FPI2,FETA2(3),FETAP2(3),
319
320 7 FRHO2,FPHI2(3),FOMEG2(3)
321
322 INTEGER MNAME(3,3,2),N4(3),I,J,ID2,ID4,I2,I4,M1,M2,IHAD1,IHAD2
323
324 LOGICAL SPIN0,SPIN1
325
326 EXTERNAL HWR,HWRUNI,HWUALF
327
328 SAVE FPI2,FETA2,FETAP2,FRHO2,FPHI2,FOMEG2,HCS,REDS,FACT,DELT,
329
330 & C1STU,C3STU
331
332 PARAMETER (EPS=1.D-20)
333
334 DATA MNAME/21,38,42,30,21,34,50,46,0,23,39,43,31,23,35,51,47,0/
335
336 DATA N4,SPIN0,SPIN1/3,3,2,.TRUE.,.TRUE./
337
338 DATA C1WVFN,FPI,FETA8,FETA1,FRHO,FPHI8,FPHI1/1.,3*0.093,3*0.107/
339
340 IF (FSTWGT) THEN
341
342 FPI2=FPI**2
343
344 CMIX=COS(ETAMIX*PIFAC/180.D0)
345
346 SMIX=SIN(ETAMIX*PIFAC/180.D0)
347
348 FETA2(1) =(+CMIX*FETA8/SQRT(TWO)-SMIX*FETA1)**2/THREE
349
350 FETA2(2) =FETA2(1)
351
352 FETA2(3) =(-CMIX*FETA8*SQRT(TWO)-SMIX*FETA1)**2/THREE
353
354 FETAP2(1)=(+SMIX*FETA8/SQRT(TWO)+CMIX*FETA1)**2/THREE
355
356 FETAP2(2)=FETAP2(1)
357
358 FETAP2(3)=(-SMIX*FETA8*SQRT(TWO)+CMIX*FETA1)**2/THREE
359
360 FRHO2=FRHO**2
361
362 CMIX=COS(PHIMIX*PIFAC/180.D0)
363
364 SMIX=SIN(PHIMIX*PIFAC/180.D0)
365
366 FPHI2(1) =(+CMIX*FPHI8/SQRT(TWO)-SMIX*FPHI1)**2/THREE
367
368 FPHI2(2) =FPHI2(1)
369
370 FPHI2(3) =(-CMIX*FPHI8*SQRT(TWO)-SMIX*FPHI1)**2/THREE
371
372 FOMEG2(1)=(+SMIX*FPHI8/SQRT(TWO)+CMIX*FPHI1)**2/THREE
373
374 FOMEG2(2)=FOMEG2(1)
375
376 FOMEG2(3)=(-SMIX*FPHI8*SQRT(TWO)+CMIX*FPHI1)**2/THREE
377
378 ENDIF
379
380 SPIN0=.NOT.(MOD(IPROC/10,10).EQ.2)
381
382 SPIN1=.NOT.(MOD(IPROC/10,10).EQ.1)
383
384 IF (GENEV) THEN
385
386 RCS=HCS*HWR()
387
388 ELSE
389
390 EVWGT=ZERO
391
392 IHAD1=1
393
394 IF (JDAHEP(1,IHAD1).NE.0) IHAD1=JDAHEP(1,IHAD1)
395
396 IHAD2=2
397
398 IF (JDAHEP(1,IHAD2).NE.0) IHAD2=JDAHEP(1,IHAD2)
399
400 PP1=PHEP(4,IHAD1)+ABS(PHEP(3,IHAD1))
401
402 PP2=PHEP(4,IHAD2)+ABS(PHEP(3,IHAD2))
403
404 XX(1)=ONE
405
406 CALL HWRPOW(ET,EJ)
407
408 EXY=EXP(HWRUNI(1,YJMIN,YJMAX))
409
410 EXY2=TWO*PP1/ET-EXY
411
412 IF (EXY2.LE.EXP(YJMIN).OR.EXY2.GE.EXP(YJMAX)) RETURN
413
414 XX(2)=PP1/(PP2*EXY*EXY2)
415
416 IF (XX(2).LE.ZERO.OR.XX(2).GE.ONE) RETURN
417
418 S=XX(2)*PP1*PP2
419
420 REDS=SQRT(S-ET*SQRT(S))
421
422 T=-HALF*PP1*ET/EXY
423
424 U=-S-T
425
426 COSTH=(T-U)/S
427
428C Set EMSCA to hard process scale (Approx ET-jet)
429
430 EMSCA=SQRT(TWO*S*T*U/(S*S+T*T+U*U))
431
432 FACT=-GEV2NB*ET*EJ*(YJMAX-YJMIN)*ALPHEM*CFFAC
433
434 & *(HWUALF(1,EMSCA)*PIFAC*C1WVFN)**2*32.D0/(THREE*S*T)
435
436 CALL HWSFUN(XX(2),EMSCA,IDHW(IHAD2),NSTRU,DISF(1,2),2)
437
438 DO 10 I=1,3
439
440 DO 10 J=1,3
441
442 10 DELT(I,J)=(QFCH(I)*U+QFCH(J)*S)**2
443
444 C1STU=-(S**2+U**2)/(T*S**2*U**2)
445
446 C3STU=-8.D0*T/(S**2*U**2)
447
448 ENDIF
449
450 HCS=ZERO
451
452 DO 50 I2=1,3
453
454C Quark initiated processes
455
456 ID2=I2
457
458 IF (DISF(ID2,2).LT.EPS) GOTO 30
459
460 DO 20 ID4=1,N4(I2)
461
462 M1=MNAME(ID2,ID4,1)
463
464 FACTR=FACT*DELT(ID2,ID4)*DISF(ID2,2)
465
466 IF (ID2.EQ.ID4) FACTR=HALF*FACTR
467
468 IF (SPIN0.AND.REDS.GT.RMASS(M1)) THEN
469
470C photon+q --> meson_0+q'
471
472 HCS=HCS+HALF*FACTR*C1STU*FPI2
473
474 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M1,ID4,1432,71,*99)
475
476 ENDIF
477
478 M2=MNAME(ID2,ID4,2)
479
480 IF (SPIN1.AND.REDS.GT.RMASS(M2)) THEN
481
482C photon+q --> meson_L+q'
483
484 HCS=HCS+FACTR*C1STU*FRHO2
485
486 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M2,ID4,1432,72,*99)
487
488C photon+q --> meson_T+q'
489
490 HCS=HCS+FACTR*C3STU*FRHO2
491
492 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M2,ID4,1432,73,*99)
493
494 ENDIF
495
496 20 CONTINUE
497
498 FACTR=FACT*DELT(I2,I2)*DISF(ID2,2)
499
500 IF (SPIN0.AND.REDS.GT.RMASS(22)) THEN
501
502C photon+q -->eta+q
503
504 HCS=HCS+HALF*FACTR*C1STU*FETA2(I2)
505
506 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(22,ID2,1432,71,*99)
507
508 ENDIF
509
510 IF (SPIN0.AND.REDS.GT.RMASS(25)) THEN
511
512C photon+q -->eta'+q
513
514 HCS=HCS+HALF*FACTR*C1STU*FETAP2(I2)
515
516 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(25,ID2,1432,71,*99)
517
518 ENDIF
519
520 IF (SPIN1.AND.REDS.GT.RMASS(56)) THEN
521
522C photon+q -->phi_L+q
523
524 HCS=HCS+FACTR*C1STU*FPHI2(I2)
525
526 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(56,ID2,1432,72,*99)
527
528C photon+q -->phi_T+q
529
530 HCS=HCS+FACTR*C3STU*FPHI2(I2)
531
532 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(56,ID2,1432,73,*99)
533
534 ENDIF
535
536 IF (SPIN1.AND.REDS.GT.RMASS(24)) THEN
537
538C photon+q -->omega_L+q
539
540 HCS=HCS+FACTR*C1STU*FOMEG2(I2)
541
542 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(24,ID2,1432,72,*99)
543
544C photon+q -->omega_T+q
545
546 HCS=HCS+FACTR*C3STU*FOMEG2(I2)
547
548 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(24,ID2,1432,73,*99)
549
550 ENDIF
551
552C Anti-quark initiated processes
553
554 30 ID2=I2+6
555
556 IF (DISF(ID2,2).LT.EPS) GOTO 50
557
558 DO 40 I4=1,N4(I2)
559
560 ID4=I4+6
561
562 FACTR=FACT*DELT(I2,I4)*DISF(ID2,2)
563
564 IF (ID2.EQ.ID4) FACTR=HALF*FACTR
565
566 M1=MNAME(I4,I2,1)
567
568 IF (SPIN0.AND.REDS.GT.RMASS(M1)) THEN
569
570C photon+qbar --> meson_0+qbar'
571
572 HCS=HCS+HALF*FACTR*C1STU*FPI2
573
574 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M1,ID4,1432,74,*99)
575
576 ENDIF
577
578 M2=MNAME(I4,I2,2)
579
580 IF (SPIN1.AND.REDS.GT.RMASS(M2)) THEN
581
582C photon+qbar --> meson_L+qbar'
583
584 HCS=HCS+FACTR*C1STU*FRHO2
585
586 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M2,ID4,1432,75,*99)
587
588C photon+qbar --> meson_T+qbar'
589
590 HCS=HCS+FACTR*C3STU*FRHO2
591
592 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(M2,ID4,1432,76,*99)
593
594 ENDIF
595
596 40 CONTINUE
597
598 FACTR=FACT*DELT(I2,I2)*DISF(ID2,2)
599
600 IF (SPIN0.AND.REDS.GT.RMASS(22)) THEN
601
602C photon+qbar -->eta+qbar
603
604 HCS=HCS+HALF*FACTR*C1STU*FETA2(I2)
605
606 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(22,ID2,1432,74,*99)
607
608 ENDIF
609
610 IF (SPIN0.AND.REDS.GT.RMASS(25)) THEN
611
612C photon+qbar -->eta'+qbar
613
614 HCS=HCS+HALF*FACTR*C1STU*FETAP2(I2)
615
616 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(25,ID2,1432,74,*99)
617
618 ENDIF
619
620 IF (SPIN1.AND.REDS.GT.RMASS(56)) THEN
621
622C photon+qbar -->phi_L+qbar
623
624 HCS=HCS+FACTR*C1STU*FPHI2(I2)
625
626 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(56,ID2,1432,75,*99)
627
628C photon+qbar -->phi_T+qbar
629
630 HCS=HCS+FACTR*C3STU*FPHI2(I2)
631
632 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(56,ID2,1432,76,*99)
633
634 ENDIF
635
636 IF (SPIN1.AND.REDS.GT.RMASS(24)) THEN
637
638C photon+qbar -->omega_L+qbar
639
640 HCS=HCS+FACTR*C1STU*FOMEG2(I2)
641
642 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(24,ID2,1432,75,*99)
643
644C photon+qbar -->omega_T+qbar
645
646 HCS=HCS+FACTR*C3STU*FOMEG2(I2)
647
648 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(24,ID2,1432,76,*99)
649
650 ENDIF
651
652 50 CONTINUE
653
654 EVWGT=HCS
655
656 RETURN
657
658C Generate event
659
660 99 IDN(1)=59
661
662 IDN(2)=ID2
663
664 IDCMF=15
665
666 CALL HWETWO
667
668C Set polarization vector
669
670 IF (IHPRO.EQ.72.OR.IHPRO.EQ.75) THEN
671
672 RHOHEP(2,NHEP-1)=ONE
673
674 ELSEIF (IHPRO.EQ.73.OR.IHPRO.EQ.76) THEN
675
676 RHOHEP(1,NHEP-1)=HALF
677
678 RHOHEP(3,NHEP-1)=HALF
679
680 ENDIF
681
682 999 END
683
684CDECK ID>, HWHPPT.
685
686*CMZ :- -12/01/93 10.12.43 by Bryan Webber
687
688*-- Author : Ian Knowles
689
690C-----------------------------------------------------------------------
691
692 SUBROUTINE HWHPPT
693
694C-----------------------------------------------------------------------
695
696C point-like photon/QCD di-jet production: mean EVWGT = sigma in nb
697
698C-----------------------------------------------------------------------
699
700 INCLUDE 'HERWIG61.INC'
701
702 DOUBLE PRECISION HWR,HWRUNI,HWUALF,EPS,RCS,PP1,PP2,ET,EJ,
703
704 & EXY,EXY2,FACTR,RS,S,T,U,CSTU,CTSU,HCS
705
706 INTEGER ID1,ID2,ID3,ID4,IHAD1,IHAD2
707
708 EXTERNAL HWR,HWRUNI,HWUALF
709
710 SAVE CSTU,CTSU,HCS,FACTR,RS
711
712 PARAMETER (EPS=1.E-9)
713
714 IHAD1=1
715
716 IF (JDAHEP(1,IHAD1).NE.0) IHAD1=JDAHEP(1,IHAD1)
717
718 IHAD2=2
719
720 IF (JDAHEP(1,IHAD2).NE.0) IHAD2=JDAHEP(1,IHAD2)
721
722 IF (GENEV) THEN
723
724 RCS=HCS*HWR()
725
726 ELSE
727
728 EVWGT=0.
729
730 PP1=PHEP(4,IHAD1)+ABS(PHEP(3,IHAD1))
731
732 PP2=PHEP(4,IHAD2)+ABS(PHEP(3,IHAD2))
733
734 XX(1)=1.
735
736 CALL HWRPOW(ET,EJ)
737
738 EXY=EXP(HWRUNI(1,YJMIN,YJMAX))
739
740 EXY2=2.*PP1/ET-EXY
741
742 IF (EXY2.LE.EXP(YJMIN).OR.EXY2.GE.EXP(YJMAX)) RETURN
743
744 XX(2)=PP1/(PP2*EXY*EXY2)
745
746 IF (XX(2).LE.ZERO.OR.XX(2).GE.ONE) RETURN
747
748 S=XX(2)*PP1*PP2
749
750 RS=.5*SQRT(S)
751
752 T=-PP1*0.5*ET/EXY
753
754 U=-S-T
755
756 COSTH=(T-U)/S
757
758C Set EMSCA to hard process scale (Approx ET-jet)
759
760 EMSCA=SQRT(2.*S*T*U/(S*S+T*T+U*U))
761
762 FACTR=-GEV2NB*0.5*EJ*(YJMAX-YJMIN)*ET*PIFAC*ALPHEM
763
764 & *HWUALF(1,EMSCA)/(S*T)
765
766 CALL HWSFUN(XX(2),EMSCA,IDHW(IHAD2),NSTRU,DISF(1,2),2)
767
768 CSTU=U/T+T/U
769
770 CTSU=-2.*CFFAC*(U/S+S/U)
771
772 ENDIF
773
774 HCS=0.
775
776 ID1=59
777
778 DO 20 ID2=1,13
779
780 IF (DISF(ID2,2).LT.EPS) GOTO 20
781
782 IF (ID2.LT.7) THEN
783
784C photon+q ---> g+q
785
786 HCS=HCS+CTSU*DISF(ID2,2)*QFCH(ID2)**2
787
788 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13,ID2,1423,51,*99)
789
790 ELSEIF (ID2.LT.13) THEN
791
792C photon+qbar ---> g+qbar
793
794 HCS=HCS+CTSU*DISF(ID2,2)*QFCH(ID2-6)**2
795
796 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13,ID2,1342,52,*99)
797
798 ELSE
799
800C photon+g ---> q+qbar
801
802 DO 10 ID3=1,6
803
804 IF (RS.GT.RMASS(ID3)) THEN
805
806 ID4=ID3+6
807
808 HCS=HCS+CSTU*DISF(ID2,2)*QFCH(ID3)**2
809
810 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID3,ID4,1423,53,*99)
811
812 ENDIF
813
814 10 CONTINUE
815
816 ENDIF
817
818 20 CONTINUE
819
820 EVWGT=FACTR*HCS
821
822 RETURN
823
824C Generate event
825
826 99 IDN(1)=ID1
827
828 IDN(2)=ID2
829
830 IDCMF=15
831
832 CALL HWETWO
833
834 999 END
835
836CDECK ID>, HWHPQS.
837
838*CMZ :- -27/03/95 13.27.22 by Mike Seymour
839
840*-- Author : Ian Knowles
841
842C-----------------------------------------------------------------------
843
844 SUBROUTINE HWHPQS
845
846C-----------------------------------------------------------------------
847
848C Compton scattering of point-like photon and (anti)quark
849
850C mean EVWGT = sigma in nb
851
852C-----------------------------------------------------------------------
853
854 INCLUDE 'HERWIG61.INC'
855
856 DOUBLE PRECISION HWR,HWRUNI,EPS,RCS,PP1,PP2,ET,EJ,EXY,EXY2,
857
858 & FACTR,S,T,U,CTSU,HCS
859
860 INTEGER ID1,ID2,IHAD1,IHAD2
861
862 EXTERNAL HWR,HWRUNI
863
864 SAVE CTSU,HCS,FACTR
865
866 PARAMETER (EPS=1.E-9)
867
868 IHAD1=1
869
870 IF (JDAHEP(1,IHAD1).NE.0) IHAD1=JDAHEP(1,IHAD1)
871
872 IHAD2=2
873
874 IF (JDAHEP(1,IHAD2).NE.0) IHAD2=JDAHEP(1,IHAD2)
875
876 IF (GENEV) THEN
877
878 RCS=HCS*HWR()
879
880 ELSE
881
882 EVWGT=0.
883
884 PP1=PHEP(4,IHAD1)+ABS(PHEP(3,IHAD1))
885
886 PP2=PHEP(4,IHAD2)+ABS(PHEP(3,IHAD2))
887
888 XX(1)=1.
889
890 CALL HWRPOW(ET,EJ)
891
892 EXY=EXP(HWRUNI(1,YJMIN,YJMAX))
893
894 EXY2=2.*PP1/ET-EXY
895
896 IF (EXY2.LE.EXP(YJMIN).OR.EXY2.GE.EXP(YJMAX)) RETURN
897
898 XX(2)=PP1/(PP2*EXY*EXY2)
899
900 IF (XX(2).LE.ZERO.OR.XX(2).GE.ONE) RETURN
901
902 S=XX(2)*PP1*PP2
903
904 T=-PP1*0.5*ET/EXY
905
906 U=-S-T
907
908 COSTH=(T-U)/S
909
910C Set EMSCA to hard process scale (Approx ET-jet)
911
912 EMSCA=SQRT(2.*S*T*U/(S*S+T*T+U*U))
913
914 FACTR=-GEV2NB*0.5*EJ*(YJMAX-YJMIN)*ET*PIFAC*ALPHEM**2/(S*T)
915
916 CALL HWSFUN(XX(2),EMSCA,IDHW(IHAD2),NSTRU,DISF(1,2),2)
917
918 CTSU=-2.*(U/S+S/U)
919
920 ENDIF
921
922 HCS=0.
923
924 ID1=59
925
926 DO 20 ID2=1,12
927
928 IF (DISF(ID2,2).LT.EPS) GOTO 20
929
930 IF (ID2.LT.7) THEN
931
932C photon+q ---> photon+q
933
934 HCS=HCS+CTSU*DISF(ID2,2)*QFCH(ID2)**4
935
936 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 59,ID2,1432,66,*99)
937
938 ELSE
939
940C photon+qbar ---> photon+qbar
941
942 HCS=HCS+CTSU*DISF(ID2,2)*QFCH(ID2-6)**4
943
944 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 59,ID2,1432,67,*99)
945
946 ENDIF
947
948 20 CONTINUE
949
950 EVWGT=FACTR*HCS
951
952 RETURN
953
954C Generate event
955
956 99 IDN(1)=ID1
957
958 IDN(2)=ID2
959
960 IDCMF=15
961
962 CALL HWETWO
963
964 999 END
965
966CDECK ID>, HWHQCD.
967
968*CMZ :- -20/05/99 12.39.45 by Kosuke Odagiri
969
970*-- Author : Bryan Webber
971
972C-----------------------------------------------------------------------
973
974 SUBROUTINE HWHQCD
975
976C-----------------------------------------------------------------------
977
978C QCD HARD 2->2 PROCESSES: MEAN EVWGT = SIGMA IN NB
979
980C-----------------------------------------------------------------------
981
982 INCLUDE 'HERWIG61.INC'
983
984 DOUBLE PRECISION HWR,HWRUNI,HWUALF,RS,EPS,HF,RCS,Z1,Z2,ET,EJ,
985
986 & FACTR,S,T,U,ST,TU,US,STU,TUS,UST,EN,RN,GFLA,AF,ASTU,ASUT,AUST,
987
988 & BF,BSTU,BSUT,BUST,BUTS,CF,CSTU,CSUT,CTSU,CTUS,DF,DSTU,DTSU,DUTS,
989
990 & DIST,HCS,UT,SU,GT,KK,KK2,YJ1INF,YJ1SUP,YJ2INF,YJ2SUP
991
992 INTEGER ID1,ID2,I
993
994 EXTERNAL HWR,HWRUNI,HWUALF
995
996 SAVE HCS,ASTU,AUST,BSTU,BSUT,BUST,BUTS,CSTU,CSUT,CTSU,CTUS,
997
998 & DSTU,DTSU,DUTS,GFLA,RCS,S,T,TU,U,US
999
1000 PARAMETER (EPS=1.E-9,HF=0.5)
1001
1002 IF (GENEV) THEN
1003
1004 RCS=HCS*HWR()
1005
1006 ELSE
1007
1008 EVWGT=0.
1009
1010 CALL HWRPOW(ET,EJ)
1011
1012 KK = ET/PHEP(5,3)
1013
1014 KK2=KK**2
1015
1016 IF (KK.GE.ONE) RETURN
1017
1018 YJ1INF = MAX( YJMIN, LOG((ONE-SQRT(ONE-KK2))/KK) )
1019
1020 YJ1SUP = MIN( YJMAX, LOG((ONE+SQRT(ONE-KK2))/KK) )
1021
1022 IF (YJ1INF.GE.YJ1SUP) RETURN
1023
1024 Z1=EXP(HWRUNI(1,YJ1INF,YJ1SUP))
1025
1026 YJ2INF = MAX( YJMIN, -LOG(TWO/KK-ONE/Z1) )
1027
1028 YJ2SUP = MIN( YJMAX, LOG(TWO/KK-Z1) )
1029
1030 IF (YJ2INF.GE.YJ2SUP) RETURN
1031
1032 Z2=EXP(HWRUNI(2,YJ2INF,YJ2SUP))
1033
1034 XX(1)=.5*(Z1+Z2)*KK
1035
1036 IF (XX(1).GE.ONE) RETURN
1037
1038 XX(2)=XX(1)/(Z1*Z2)
1039
1040 IF (XX(2).GE.ONE) RETURN
1041
1042 COSTH=(Z1-Z2)/(Z1+Z2)
1043
1044 S=XX(1)*XX(2)*PHEP(5,3)**2
1045
1046 RS=HF*SQRT(S)
1047
1048 DO 3 I=1,NFLAV
1049
1050 IF (RS.LT.RMASS(I)) GOTO 4
1051
1052 3 CONTINUE
1053
1054 I=NFLAV+1
1055
1056 4 MAXFL=I-1
1057
1058 IF (MAXFL.EQ.0) CALL HWWARN('HWHQCD',100,*999)
1059
1060C
1061
1062 T=-HF*S*(1.-COSTH)
1063
1064 U=-S-T
1065
1066C---SET EMSCA TO HARD PROCESS SCALE (APPROX ET-JET)
1067
1068 EMSCA=SQRT(2.*S*T*U/(S*S+T*T+U*U))
1069
1070 FACTR = GEV2NB*.5*PIFAC*EJ*ET*(HWUALF(1,EMSCA)/S)**2
1071
1072 & * (YJ1SUP-YJ1INF)*(YJ2SUP-YJ2INF)
1073
1074 CALL HWSGEN(.FALSE.)
1075
1076C
1077
1078 ST=S/T
1079
1080 TU=T/U
1081
1082 US=U/S
1083
1084 STU=TU/US
1085
1086 TUS=US/ST
1087
1088 UST=ST/TU
1089
1090C
1091
1092 EN=CAFAC
1093
1094 RN=CFFAC/EN
1095
1096 GFLA=HF*FLOAT(MAXFL)/(EN*RN)**2
1097
1098 AF=FACTR*RN
1099
1100 ASTU=AF*(1.-2.*UST)
1101
1102 ASUT=AF*(1.-2.*STU)
1103
1104 AUST=AF*(1.-2.*TUS)
1105
1106C-----------------------------------------------------------------------
1107
1108C---Colour decomposition modifications below (KO)
1109
1110C-----------------------------------------------------------------------
1111
1112 BF=HF-AF/EN/TUS/(ASTU+ASUT)
1113
1114 BSTU=BF*ASTU
1115
1116 BSUT=BF*ASUT
1117
1118 BF=ONE-TWO*AF/EN/STU/(AUST+ASTU)
1119
1120 BUST=BF*AUST
1121
1122 BUTS=BF*ASTU
1123
1124C-----------------------------------------------------------------------
1125
1126C BF=2.*AF/EN
1127
1128C BSTU=HF*(ASTU+BF*ST)
1129
1130C BSUT=HF*(ASUT+BF/US)
1131
1132C BUST=AUST+BF*US
1133
1134C BUTS=ASTU+BF/TU
1135
1136C-----------------------------------------------------------------------
1137
1138 CF=AF*EN
1139
1140 CSTU=(CF*(RN-TUS))/TU
1141
1142 CSUT=(CF*(RN-TUS))*TU
1143
1144 CTSU=(FACTR*(UST-RN))*US
1145
1146 CTUS=(FACTR*(UST-RN))/US
1147
1148 DF=HF*FACTR/RN
1149
1150 DSTU=DF*(1.+1./TUS-STU-UST)
1151
1152 DTSU=DF*(1.+1./UST-STU-TUS)
1153
1154 DUTS=DF*(1.+1./STU-UST-TUS)
1155
1156 ENDIF
1157
1158C
1159
1160 HCS=0.
1161
1162 DO 6 ID1=1,13
1163
1164 IF (DISF(ID1,1).LT.EPS) GOTO 6
1165
1166 DO 5 ID2=1,13
1167
1168 IF (DISF(ID2,2).LT.EPS) GOTO 5
1169
1170 DIST=DISF(ID1,1)*DISF(ID2,2)
1171
1172 IF (ID1.LT.7) THEN
1173
1174C---QUARK FIRST
1175
1176 IF (ID2.LT.7) THEN
1177
1178 IF (ID1.NE.ID2) THEN
1179
1180 HCS=HCS+ASTU*DIST
1181
1182 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421, 3,*9)
1183
1184 ELSE
1185
1186 HCS=HCS+BSTU*DIST
1187
1188 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421, 1,*9)
1189
1190 HCS=HCS+BSUT*DIST
1191
1192 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,4312, 2,*9)
1193
1194 ENDIF
1195
1196 ELSEIF (ID2.NE.13) THEN
1197
1198 IF (ID2.NE.ID1+6) THEN
1199
1200 HCS=HCS+ASTU*DIST
1201
1202 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3142, 9,*9)
1203
1204 ELSE
1205
1206 HCS=HCS+FLOAT(MAXFL-1)*AUST*DIST
1207
1208 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(-ID1, 0,2413, 4,*9)
1209
1210 HCS=HCS+BUTS*DIST
1211
1212 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3142, 5,*9)
1213
1214 HCS=HCS+BUST*DIST
1215
1216 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413, 6,*9)
1217
1218 HCS=HCS+CSTU*DIST
1219
1220 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13, 13,2413, 7,*9)
1221
1222 HCS=HCS+CSUT*DIST
1223
1224 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13, 13,2341, 8,*9)
1225
1226 ENDIF
1227
1228 ELSE
1229
1230 HCS=HCS+CTSU*DIST
1231
1232 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3142,10,*9)
1233
1234 HCS=HCS+CTUS*DIST
1235
1236 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421,11,*9)
1237
1238 ENDIF
1239
1240 ELSEIF (ID1.NE.13) THEN
1241
1242C---QBAR FIRST
1243
1244 IF (ID2.LT.7) THEN
1245
1246 IF (ID1.NE.ID2+6) THEN
1247
1248 HCS=HCS+ASTU*DIST
1249
1250 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413,17,*9)
1251
1252 ELSE
1253
1254 HCS=HCS+FLOAT(MAXFL-1)*AUST*DIST
1255
1256 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(-ID1, 0,3142,12,*9)
1257
1258 HCS=HCS+BUTS*DIST
1259
1260 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413,13,*9)
1261
1262 HCS=HCS+BUST*DIST
1263
1264 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3142,14,*9)
1265
1266 HCS=HCS+CSTU*DIST
1267
1268 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13, 13,3142,15,*9)
1269
1270 HCS=HCS+CSUT*DIST
1271
1272 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 13, 13,4123,16,*9)
1273
1274 ENDIF
1275
1276 ELSEIF (ID2.NE.13) THEN
1277
1278 IF (ID1.NE.ID2) THEN
1279
1280 HCS=HCS+ASTU*DIST
1281
1282 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,4312,20,*9)
1283
1284 ELSE
1285
1286 HCS=HCS+BSTU*DIST
1287
1288 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,4312,18,*9)
1289
1290 HCS=HCS+BSUT*DIST
1291
1292 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421,19,*9)
1293
1294 ENDIF
1295
1296 ELSE
1297
1298 HCS=HCS+CTSU*DIST
1299
1300 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413,21,*9)
1301
1302 HCS=HCS+CTUS*DIST
1303
1304 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,4312,22,*9)
1305
1306 ENDIF
1307
1308 ELSE
1309
1310C---GLUON FIRST
1311
1312 IF (ID2.LT.7) THEN
1313
1314 HCS=HCS+CTSU*DIST
1315
1316 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413,23,*9)
1317
1318 HCS=HCS+CTUS*DIST
1319
1320 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421,24,*9)
1321
1322 ELSEIF (ID2.LT.13) THEN
1323
1324 HCS=HCS+CTSU*DIST
1325
1326 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3142,25,*9)
1327
1328 HCS=HCS+CTUS*DIST
1329
1330 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,4312,26,*9)
1331
1332 ELSE
1333
1334 HCS=HCS+GFLA*CSTU*DIST
1335
1336 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 0, 0,2413,27,*9)
1337
1338 HCS=HCS+GFLA*CSUT*DIST
1339
1340 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP( 0, 0,4123,28,*9)
1341
1342 HCS=HCS+DTSU*DIST
1343
1344 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2341,29,*9)
1345
1346 HCS=HCS+DSTU*DIST
1347
1348 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,3421,30,*9)
1349
1350 HCS=HCS+DUTS*DIST
1351
1352 IF (GENEV.AND.HCS.GT.RCS) CALL HWHQCP(ID1,ID2,2413,31,*9)
1353
1354 ENDIF
1355
1356 ENDIF
1357
1358 5 CONTINUE
1359
1360 6 CONTINUE
1361
1362 EVWGT=HCS
1363
1364 RETURN
1365
1366C---GENERATE EVENT
1367
1368 9 IDN(1)=ID1
1369
1370 IDN(2)=ID2
1371
1372 IDCMF=15
1373
1374 CALL HWETWO
1375
1376 IF (AZSPIN) THEN
1377
1378C Calculate coefficients for constructing spin density matrices
1379
1380 IF (IHPRO.EQ.7 .OR.IHPRO.EQ.8 .OR.
1381
1382 & IHPRO.EQ.15.OR.IHPRO.EQ.16) THEN
1383
1384C qqbar-->gg or qbarq-->gg
1385
1386 UT=1./TU
1387
1388 GCOEF(1)=UT+TU
1389
1390 GCOEF(2)=-2.
1391
1392 GCOEF(3)=0.
1393
1394 GCOEF(4)=0.
1395
1396 GCOEF(5)=GCOEF(1)
1397
1398 GCOEF(6)=UT-TU
1399
1400 GCOEF(7)=-GCOEF(6)
1401
1402 ELSEIF (IHPRO.EQ.10.OR.IHPRO.EQ.11.OR.
1403
1404 & IHPRO.EQ.21.OR.IHPRO.EQ.22.OR.
1405
1406 & IHPRO.EQ.23.OR.IHPRO.EQ.24.OR.
1407
1408 & IHPRO.EQ.25.OR.IHPRO.EQ.26) THEN
1409
1410C qg-->qg or qbarg-->qbarg or gq-->gq or gqbar-->gqbar
1411
1412 SU=1./US
1413
1414 GCOEF(1)=-(SU+US)
1415
1416 GCOEF(2)=0.
1417
1418 GCOEF(3)=2.
1419
1420 GCOEF(4)=0.
1421
1422 GCOEF(5)=SU-US
1423
1424 GCOEF(6)=GCOEF(1)
1425
1426 GCOEF(7)=-GCOEF(5)
1427
1428 ELSEIF (IHPRO.EQ.27.OR.IHPRO.EQ.28) THEN
1429
1430C gg-->qqbar
1431
1432 UT=1./TU
1433
1434 GCOEF(1)=TU+UT
1435
1436 GCOEF(2)=-2.
1437
1438 GCOEF(3)=0.
1439
1440 GCOEF(4)=0.
1441
1442 GCOEF(5)=GCOEF(1)
1443
1444 GCOEF(6)=TU-UT
1445
1446 GCOEF(7)=-GCOEF(6)
1447
1448 ELSEIF (IHPRO.EQ.29.OR.IHPRO.EQ.30.OR.
1449
1450 & IHPRO.EQ.31) THEN
1451
1452C gg-->gg
1453
1454 GT=S*S+T*T+U*U
1455
1456 GCOEF(2)=2.*U*U*T*T
1457
1458 GCOEF(3)=2.*S*S*U*U
1459
1460 GCOEF(4)=2.*S*S*T*T
1461
1462 GCOEF(1)=GT*GT-GCOEF(2)-GCOEF(3)-GCOEF(4)
1463
1464 GCOEF(5)=GT*(GT-2.*S*S)-GCOEF(2)
1465
1466 GCOEF(6)=GT*(GT-2.*T*T)-GCOEF(3)
1467
1468 GCOEF(7)=GT*(GT-2.*U*U)-GCOEF(4)
1469
1470 ELSE
1471
1472 CALL HWVZRO(7,GCOEF)
1473
1474 ENDIF
1475
1476 ENDIF
1477
1478 999 END
1479
1480CDECK ID>, HWHQCP.
1481
1482*CMZ :- -26/04/91 10.18.57 by Bryan Webber
1483
1484*-- Author : Bryan Webber
1485
1486C-----------------------------------------------------------------------
1487
1488 SUBROUTINE HWHQCP(ID3,ID4,IPERM,IHPR,*)
1489
1490C-----------------------------------------------------------------------
1491
1492C IDENTIFIES HARD SUBPROCESS
1493
1494C-----------------------------------------------------------------------
1495
1496 INCLUDE 'HERWIG61.INC'
1497
1498 INTEGER HWRINT,ID3,ID4,IPERM,IHPR,ND3
1499
1500 EXTERNAL HWRINT
1501
1502 IHPRO=IHPR
1503
1504 IF (ID3.GT.0) THEN
1505
1506 IDN(3)=ID3
1507
1508 IDN(4)=ID4
1509
1510 ELSE
1511
1512 ND3=-ID3
1513
1514 IF (ID3.GT.-7) THEN
1515
1516 1 IDN(3)=HWRINT(1,MAXFL)
1517
1518 IF (IDN(3).EQ.ND3) GOTO 1
1519
1520 IDN(4)=IDN(3)+6
1521
1522 ELSE
1523
1524 2 IDN(3)=HWRINT(1,MAXFL)+6
1525
1526 IF (IDN(3).EQ.ND3) GOTO 2
1527
1528 IDN(4)=IDN(3)-6
1529
1530 ENDIF
1531
1532 ENDIF
1533
1534 ICO(1)=IPERM/1000
1535
1536 ICO(2)=IPERM/100-10*ICO(1)
1537
1538 ICO(3)=IPERM/10 -10*(IPERM/100)
1539
1540 ICO(4)=IPERM -10*(IPERM/10)
1541
1542 RETURN 1
1543
1544 END