* * $Id: seclf5.F,v 1.1.1.1 1995/10/24 10:21:59 cernlib Exp $ * * $Log: seclf5.F,v $ * Revision 1.1.1.1 1995/10/24 10:21:59 cernlib * Geant * * #include "geant321/pilot.h" *CMZ : 3.21/04 23/02/95 14.46.01 by S.Giani *-- Author : SUBROUTINE SECLF5(FSE,IFSE,EX,U,E) C THIS ROUTINE SAMPLES AN EXIT ENERGY FROM C A GENERAL EVAPORATION SPECTRUM #include "geant321/minput.inc" DIMENSION FSE(*),IFSE(*) SAVE C DETERMINE THETA IP=1 NR=IFSE(IP) NE=IFSE(IP+1) IP=2*NR+IP EMAX=E-U R=FLTRNF(0) DO 10 I=1,NE IP=IP+2 IF(E.LE.FSE(IP))GO TO 20 10 CONTINUE GO TO 80 20 IF(I.EQ.1)GO TO 90 C DETERMINE THE INTERPOLATING SCHEME CALL INTSCH(IFSE,I,IS,NR) E2=FSE(IP) E1=FSE(IP-2) CALL INTERP(E,THETA,E1,FSE(IP-1),E2,FSE(IP+1),IS) IP=IP+2+(NE-I)*2 C DETERMINE X 30 NF=IFSE(IP+1) NR=IFSE(IP) IPR=IP IP=IP+1+2*NR PROB=0. DO 60 I=1,NF N=IP+2*I PROB1=PROB CALL INTSCH(IFSE(IPR),I,IS,NR) IF(I.EQ.1)GO TO 60 IF(IS.EQ.1)GO TO 40 IF(IS.GT.2)GO TO 100 PROB=PROB+(FSE(N)+FSE(N-2))*(FSE(N-1)-FSE(N-3))/2. GO TO 50 40 PROB=PROB+FSE(N-2)*(FSE(N-1)-FSE(N-3)) 50 CONTINUE IF(R.LE.PROB)GO TO 70 60 CONTINUE 70 X=FSE(N-3)+(R-PROB1)*(FSE(N-1)-FSE(N-3))/(PROB-PROB1) C SELECT THE EXIT ENERGY FROM THE GENERAL EVAPORATION SPECTRUM EX=THETA*X IF(EX.LE.EMAX)RETURN #if defined(CERNLIB_MDEBUG) WRITE(IOUT,10000)EX,EMAX 10000 FORMAT(' MICAP: WARNING-EX,EMAX=',1P2E13.5,' IN SECLF5') #endif EX=EMAX RETURN C INCIDENT ENERGY IS ABOVE THE LAST INCIDENT ENERGY GIVEN C USE THE LAST DISTRIBUTION 80 THETA=FSE(IP+1) IP=IP+2 GO TO 30 C INCIDENT ENERGY IS BELOW THE FIRST INCIDENT ENERGY GIVEN C USE THE FIRST DISTRIBUTION 90 THETA=FSE(IP+1) IP=IP+2*(NE-I)+2 GO TO 30 100 WRITE(IOUT,10100)IS 10100 FORMAT(' MICAP: INTERPOLATION SCHEME=',I3,' IN SECLF5') WRITE(6,*) ' CALOR: ERROR in SECLF5 =====> STOP ' STOP END