*
* $Id: genev.F,v 1.1.1.1 1996/03/22 16:42:46 mclareni Exp $
*
* $Log: genev.F,v $
* Revision 1.1.1.1  1996/03/22 16:42:46  mclareni
* Phtools
*
*
      SUBROUTINE GENEV
C--      CALLED FROM - MAIN PROGRAM
C --     PURPOSE - GENERATES EVENT (PCM AND WT ) USING AMASS AND TECM
      COMMON/UTIL/PLACE(150),WEIGHT(150),NP,TECM,AMASS(18),PCM(5,18),WT
C--       UTIL2 CONTAINS FURTHER NUMBERS USEFUL TO USER,
C--            NAMELY EFF.MASSOF PARTICLES 1-I, M**2, AND DELTA SQ.
C--       M**2 AND DELTA SQ. ARE INDEFINITE IF RAUBOLD-LYNCH METHOD IS U
      COMMON/UTIL2/ EMM(18), EMMS(18), D2(18)
     +/BEAM  / PRM       ,TARGM    ,PBM      ,ETA      ,GAM      ,DPBM
     +/KGENEV/ KGENEV
     +/KAJBYK/ KAJBYK    ,A(18)    ,UCUT(18) ,BCUT(18) ,NN
     +/SHUFFL/ RNO(50)   ,NTNM4    ,NTM2     ,NTM1     ,QQCONS(50) ,TIR
C--       PCM1 IS LINEAR EQUIV. OF PCM TO AVOID DOUBLE INDICES
      DIMENSION EM(18),PD(18),EMS(18),SM(18),FFQ(18),PCM1(90)
      EQUIVALENCE (NT,NP),(AMASS(1),EM(1)),(PCM1(1),PCM(1,1))
C FFQ(N)=PI * (TWOPI)**(N-2) / (N-2)FACTORIAL
      DATA  FFQ/0.,3.141592, 19.73921, 62.01255, 129.8788, 204.0131,
     +                       256.3704, 268.4705, 240.9780, 189.2637,
     +                       132.1308,  83.0202,  47.4210,  24.8295,
     +                        12.0006,   5.3858,   2.2560,   0.8859/
      DATA ETC, INIT, TWOPI / 0. , 0 , 6.2831853073 /
      DDLAM(X,Y,Z) = SQRT((X-Y-Z)**2 - 4.*Z*Y)
C*IA  Y,Z
      DATA Y,Z/2*0/
      IF(TECM-ETC) 115,455,115
C        INITIALIZATION
  115 IF(INIT.GT.0) GOTO 300
      INIT=1
      NTM1=NT-1
      NTM2=NT-2
      NTP1=NT+1
      NTNM4 = 3*NT - 4
      EMM(1)=EM(1)
      TM=0.0
      TAM=PRM**2
      TARGM2 = TARGM**2
      PRMINV = 1./(2.*PRM)
      DO 200 I=1,NT
      EMS(I)=EM(I)**2
      TM=TM+EM(I)
 200  SM(I)=TM
C        CONSTANTS DEPENDING ON TECM
  300 CONTINUE
      TECMTM=TECM-TM
      ETC=TECM
      EMM(NT)=TECM
      IF(KAJBYK.LE.0) GOTO 320
      TECM2 = TECM**2
      DDL11 = DDLAM(TECM2,TARGM2,TAM)
      QX1 = DDL11*PRMINV
      ETE = QX1/TECM
      GEM=SQRT(1.+ETE**2)
  320 CONTINUE
      IF(KGENEV.GT.1) GOTO 400
C        CONSTANT CROSS SECTION AS FUNCTION OF TECM
      EMMAX=TECMTM+EM(1)
      EMMIN=0.0
      WTMAX=1.0
      DO 350 I=2,NT
      EMMIN=EMMIN+EM(I-1)
      EMMAX=EMMAX+EM(I)
  350 WTMAX=WTMAX*PDK(EMMAX,EMMIN,EM(I))
      WTMAXQ=1.0/WTMAX
      GOTO 455
C--      FERMI ENERGY DEPENDENCE FOR CROSS SECTION
  400 WTMAXQ=TECMTM**NTM2*FFQ(NT) / TECM
C        CALCULATION OF WT BASED ON EFFECTIVE MASSES EMM
C--       RANGNR FILLS RNO WITH 3*NT-4 RANDOM NUMBERS,
  455 CALL RANGNR
C--            OF WHICH THE FIRST NT-2 ARE ORDERED.
      IF(NTM2) 900,509,460
  460 CONTINUE
      DO 508 J=2,NTM1
  508 EMM(J)=RNO(J-1)*(TECMTM)+SM(J)
  509 WT=WTMAXQ
      IR=NTM2
      DO 530 I=1,NTM1
      PD(I)=PDK(EMM(I+1),EMM(I),EM(I+1))
  530 WT=WT*PD(I)
C--       KAJBYK IS FLAG (=1 IF KAJANTIE-BYCKLING METHOD TO BE USED)
      IF(KAJBYK.GT.0) GOTO 800
C--       COMPLETE SPECIFICATION OF EVENT (RAUBOLD-LYNCH METHOD)
      PCM(1,1)=0.0
      PCM(2,1)=PD(1)
      PCM(3,1)=0.0
      DO 570 I=2,NT
      PCM(1,I)=0.0
      PCM(2,I) = -PD(I-1)
      PCM(3,I)=0.0
      IR=IR+1
      BANG=TWOPI*RNO(IR)
      CB=COS(BANG)
      SB=SIN(BANG)
      IR=IR+1
      C=2.0*RNO(IR)-1.0
      S=SQRT(1.0-C*C)
      IF(I.EQ.NT) GOTO 1567
      ESYS=SQRT(PD(I)**2+EMM(I)**2)
      BETA=PD(I)/ESYS
      GAMA=ESYS/EMM(I)
      DO 568 J=1,I
      NDX = 5*J - 5
      AA= PCM1(NDX+1)**2 + PCM1(NDX+2)**2 + PCM1(NDX+3)**2
      PCM1(NDX+5) = SQRT(AA)
      PCM1(NDX+4) = SQRT(AA+EMS(J))
      CALL ROTES2(C,S,CB,SB,PCM,J)
      PSAVE = GAMA*(PCM(2,J)+BETA*PCM(4,J))
  568 PCM(2,J)=PSAVE
      GOTO 570
 1567 DO 1568 J=1,I
      AA=PCM(1,J)**2 + PCM(2,J)**2 + PCM(3,J)**2
      PCM(5,J)=SQRT(AA)
      PCM(4,J)=SQRT(AA+EMS(J))
      CALL ROTES2(C,S,CB,SB,PCM,J)
 1568 CONTINUE
  570 CONTINUE
  900 CONTINUE
      RETURN
C--       COMPLETE SPECIFICATION OF EVENT (BYCKLING-KAJANTIE METHOD)
  800 CONTINUE
      DO 805 I=1,NP
  805 EMMS(I) = EMM(I)**2
      PAM = TARGM2
      SINF1=0.
      COSF1=1.
      COST1=1.
      SINT1=0.
      QY=0.
      QZ=0.
      QX = QX1
      DDL1 = DDL11
C--       INDICES RUN BACKWARDS SINCE EVENT GENERATED IN ANTILAB
      DO 840 I=1,NTM1
      NTP1I = NTP1 - I
      EMNP1 = EMMS(NTP1I)
      EEEM = EMS(NTP1I)
      NTMI = NT - I
      EMN = EMMS(NTMI)
      PPAM=PAM
      X = (-EMNP1 +PPAM -TAM) * (EMNP1 +EMN -EEEM)
      XX = DDL1 * PD(NTMI) * EMM(NTP1I) * 2.0
      TT1 = EMN + TAM + (X-XX)/(2.0*EMNP1)
      TT2 = EMN + TAM + (X+XX)/(2.0*EMNP1)
      WP = TT2-TT1
C--       BCUT AND UCUT HAVE OPPOSITE SIGN CONVENTION OF TT1,TT2
      IF(BCUT(NTMI).NE.0.) TT2 = MIN(TT2,-BCUT(NTMI))
      IF(UCUT(NTMI).NE.0.) TT1 = MAX(TT1,-UCUT(NTMI))
      IF(TT2.LE.TT1) GOTO 455
      ANN = A(NTMI)
      IF(ANN.EQ.0.) GOTO 810
      TT1 = EXP(ANN*TT1)
      TT2 = EXP(ANN*TT2)
      TT2 = TT2 - TT1
      PAM = TT1 + RNO(IR+2)*TT2
      WT = WT / (PAM*ANN)
      WT = WT * TT2 / WP
      PAM = LOG(PAM)/ANN
      GOTO 820
  810 TT2 = TT2 - TT1
      PAM = TT1 + RNO(IR+2)*TT2
      WT = WT * TT2 / WP
  820 D2(NTP1I-1) = -PAM
      DDL2 = DDLAM(EMN,PAM,TAM)
      COST = ((TAM+PPAM-EMNP1) *(TAM+PAM-EMN) - 2.*TAM*(PPAM+PAM-EEEM))
     +      /(DDL1*DDL2)
      COST22 = COST*COST
      IF(COST22.GT.1.) GOTO 455
      SINT = SQRT(1.-COST22)
      BANG=TWOPI*RNO(IR+1)
      SINDEP=SIN(BANG)
      COSDEP=COS(BANG)
      COST2=COST1*COST+SINT1*SINT*COSDEP
      COST22 = COST2*COST2
      IF(COST22.GT.1.) GOTO 455
      SINT2= SQRT(1.-COST22)
      VK=SINT1*SINT2
      IF(VK.EQ.0.) GOTO 830
      COSDEP=(COST-COST1*COST2)/VK
      SINDEP=SINDEP*SINT/SINT2
  830 COSF2=COSF1*COSDEP-SINF1*SINDEP
      SINF2=SINF1*COSDEP+COSF1*SINDEP
      QP = DDL2*PRMINV
      QYY=QP*SINT2*SINF2
      QZZ=QP*SINT2*COSF2
      QXX=QP*COST2
      NDX = 5*NTMI
      PCM1(NDX+1) = QX-QXX
      PCM1(NDX+2) = QY-QYY
      PCM1(NDX+3) = QZ-QZZ
      DDL1=DDL2
      COST1=COST2
      SINT1=SINT2
      COSF1=COSF2
      SINF1=SINF2
      QX=QXX
      QY=QYY
      QZ=QZZ
      IR=IR+2
  840 CONTINUE
      PCM(1,1)=QX
      PCM(2,1)=QY
      PCM(3,1)=QZ
C--       TRANSFORM FROM ANTILAB TO CM SYSTEM
      DO 855 J= 1, NP
      NDX = 5*J - 5
      PXI = PCM1(NDX+1)
      AA = PXI**2 + PCM1(NDX+2)**2 + PCM1(NDX+3)**2
      EI = SQRT(AA + EMS(J))
      PCM1(NDX+1) = ETE*EI - GEM*PXI
      PCM1(NDX+4) = GEM*EI - ETE*PXI
  855 PCM1(NDX+5) = SQRT(PCM1(NDX+4)**2 - EMS(J))
      RETURN
      END
