1 #include "isajet/pilot.h"
3 C
4 C          Generate masses for uniform NADD-body phase space in DECPS2.
5 C          Auxiliary routine for DECAY.
6 C
7 #if defined(CERNLIB_IMPNONE)
8       IMPLICIT NONE
9 #endif
10 C
11 #include "isajet/itapes.inc"
12 #include "isajet/partcl.inc"
13 C
15       REAL PGEN(5,5)
16       REAL REDUCE(5),RND(5)
17       REAL RANF,PCM,DBLPCM
18       REAL WTMAX,SUM1,SUM2,SUM,RNEW,WT,A,B,C
20 C
21 C          Function definitions.
22 C
23 #if defined(CERNLIB_SINGLE)
24       PCM(A,B,C)=SQRT((A-B-C)*(A+B+C)*(A-B+C)*(A+B-C))/(2.*A)
25 #endif
26 #if defined(CERNLIB_DOUBLE)
27       PCM(A,B,C)=DBLPCM(A,B,C)
28 #endif
29 C
30       DATA REDUCE/1.,1.,2.,5.,15./
31 C
32 C          Calculate maximum phase-space weight.
33 C
37       SUM=0
38       DO 100 I=1,NADD
39         SUM=SUM+PPTCL(5,NPTCL+I)
40 100   CONTINUE
41       SUM1=PGEN(5,1)
42       SUM2=SUM-PPTCL(5,NPTCL+1)
43       DO 110 I=1,NADD1
44         WTMAX=WTMAX*PCM(SUM1,SUM2,PPTCL(5,NPTCL+I))
45         SUM1=SUM1-PPTCL(5,NPTCL+I)
46         SUM2=SUM2-PPTCL(5,NPTCL+I+1)
47 110   CONTINUE
48 C
49 C          Generate masses for uniform NADD-body phase space.
50 C
51 200   CONTINUE
52       RND(1)=1.
53       DO 210 I=2,NADD1
54         RNEW=RANF()
55         I1=I-1
56         DO 220 JJ1=1,I1
57           J=I-JJ1
58           JSAVE=J+1
59           IF(RNEW.LE.RND(J)) GO TO 210
60           RND(JSAVE)=RND(J)
61 220     CONTINUE
62 210   RND(JSAVE)=RNEW