GEANT321/neutron/secegy.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1995/10/24 10:21:58  cernlib
   6 * Geant
   7 *
   8 *
   9 #include "geant321/pilot.h"
  10 *CMZ :  3.21/04 23/02/95  14.46.01  by  S.Giani
  11 *-- Author :
  12       SUBROUTINE SECEGY(EX,FSE,E,IFSE)
  13 C       THIS ROUTINE SELECTS A PARTIAL ENERGY DISTRIBUTION
  14 C       TO SAMPLE THE EXIT ENERGY FROM
  15 #include "geant321/minput.inc"
  16       DIMENSION FSE(*),IFSE(*)
  17       SAVE
  18       EX = 0.0
  19       IPP=1
  20       N=1
  21       IP=1
  22       R=FLTRNF(0)
  23       NK=IFSE(IP)
  24       PROB=0.
  25    10 IP=IP+1
  26       LF=IFSE(IP)
  27       IP=IP+1
  28 C       TEMP FIX UP
  29       U=FSE(IP)
  30       IF(LF.EQ.11)U=FLOAT(IFSE(IP))
  31       IP=IP+1
  32       NR=IFSE(IP)
  33       IPR=IP
  34       IP=IP+1
  35       NP=IFSE(IP)
  36       IP=IP+2*NR
  37    20 CONTINUE
  38       DO 30 I=1,NP
  39          IP=IP+2
  40 C       IF E IS LESS THAN THE LOWEST ENERGY OF THE MESH, THEN THE
  41 C       PROBABILITY WILL EQUAL ZERO FOR SELECTING THAT DISTRIBUTION
  42          IF(E.LT.FSE(IP-1))GO TO 50
  43    30 CONTINUE
  44 C       TRY THE NEXT PARTIAL DISTRIBUTION
  45    40 N=N+1
  46       IF(N.GT.NK)GO TO 170
  47       IF(LF.EQ.1)GO TO 100
  48       IF(LF.EQ.5)GO TO 120
  49       IF((LF.EQ.7).OR.(LF.EQ.9))GO TO 130
  50       GO TO 140
  51    50 IF(I.NE.1)GO TO 70
  52       IF(E+CADIG(E).LT.FSE(IP-1))GO TO 60
  53       E=E+CADIG(E)
  54       IP=IP-2
  55       GO TO 20
  56    60 CONTINUE
  57       IP=IP+(NP-1)*2
  58       GO TO 40
  59 C       DETERMINE THE INTERPOLATING SCHEME
  60    70 CONTINUE
  61       DO 80 J=1,NR
  62          J1=IPR+2*J
  63          IF(I.LE.IFSE(J1))GO TO 90
  64    80 CONTINUE
  65    90 IS=IFSE(J1+1)
  66       CALL INTERP(E,P,FSE(IP-3),FSE(IP-2),FSE(IP-1),FSE(IP),IS)
  67       PROB=PROB+P
  68       IF(R.LE.PROB)GO TO 150
  69       IP=IP+2*(NP-I)
  70       GO TO 40
  71 C       SKIP THE DATA FOR LF EQUAL ONE
  72   100 IP=IP+1
  73       NR=IFSE(IP)
  74       NE=IFSE(IP+1)
  75       IP=IP+2*NR+1
  76       DO 110 I=1,NE
  77          IP=IP+2
  78          NR=IFSE(IP)
  79          IP=IP+1
  80          NP=IFSE(IP)
  81          IP=IP+2*NR+2*NP
  82   110 CONTINUE
  83       GO TO 10
  84 C       SKIP THE DATA FOR LF EQUAL FIVE
  85   120 IP=IP+1
  86       NR=IFSE(IP)
  87       NE=IFSE(IP+1)
  88       IP=IP+2*NR+1
  89       IP=IP+2*NE
  90       IP=IP+1
  91       NR=IFSE(IP)
  92       NF=IFSE(IP+1)
  93       IP=IP+2*NF+2*NR+1
  94       GO TO 10
  95 C       SKIP THE DATA FOR LF EQUAL SEVEN, AND LF EQUAL NINE
  96   130 IP=IP+1
  97       NR=IFSE(IP)
  98       NE=IFSE(IP+1)
  99       IP=IP+2*NR+1
 100       IP=IP+2*NE
 101       GO TO 10
 102 C       SKIP THE DATA FOR LF EQUAL ELEVEN
 103   140 IP=IP+1
 104       NR=IFSE(IP)
 105       NE=IFSE(IP+1)
 106       IP=IP+2*NR+1
 107       IP=IP+2*NE
 108       IP=IP+1
 109       NR=IFSE(IP)
 110       NE=IFSE(IP+1)
 111       IP=IP+2*NR+1
 112       IP=IP+2*NE
 113       GO TO 10
 114 C       NOW SELECT THE SECONDARY ENERGY FROM THE CHOSEN DISTRIBUTION
 115   150 IP=IP+2*(NP-I)
 116   160 CONTINUE
 117       IF(LF.EQ.1)CALL SECLF1(FSE(IP+1),IFSE(IP+1),EX,U,E)
 118       IF(LF.EQ.5)CALL SECLF5(FSE(IP+1),IFSE(IP+1),EX,U,E)
 119       IF(LF.EQ.7)CALL SECLF7(FSE(IP+1),IFSE(IP+1),EX,U,E)
 120       IF(LF.EQ.9)CALL SECLF9(FSE(IP+1),IFSE(IP+1),EX,U,E)
 121       IF(LF.EQ.11)CALL SECL11(FSE(IP+1),IFSE(IP+1),EX,U,E)
 122       RETURN
 123   170 CONTINUE
 124 #if defined(CERNLIB_MDEBUG)
 125       WRITE(IOUT,10000)R,PROB,E
 126 10000 FORMAT(' MICAP: WARNING-SECONDARY ENERGY DISTRIBUTION NOT ',
 127      +       'CHOSEN IN SECEGY',1P3E11.4)
 128 #endif
 129 C       TEMP CARD
 130       LF=IFSE(IPP+1)
 131       U=FSE(IPP+2)
 132       IF(LF.EQ.11)U=FLOAT(IFSE(IPP+2))
 133       NR=IFSE(IPP+3)
 134       NP=IFSE(IPP+4)
 135       IP=2*NR+2*NP+5
 136       GO TO 160
 137       END