PDF/spdf/gamfun.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.2  1996/10/30 08:28:27  cernlib
   6 * Version 7.04
   7 *
   8 * Revision 1.1.1.1  1996/04/12 15:29:22  plothow
   9 * Version 7.01
  10 *
  11 *
  12 #include "pdf/pilot.h"
  13 C-----------------------------------------------------------------------
  14 C     Gamma function : See Abramowitz, page 257, form. 6.4.40
  15 C
  16 #if defined(CERNLIB_SINGLE)
  17       FUNCTION GAMFUN(Y)
  18 #endif
  19 #if defined(CERNLIB_DOUBLE)
  20       DOUBLE PRECISION FUNCTION GAMFUN(Y)
  21 #endif
  22 #include "pdf/impdp.inc"
  23 #include "pdf/expdp.inc"
  24      +     Y,R,S,T,AFSPL,X,
  25      +     COEF(10),PI,ZEROD,HALFD,ONED,TWOD,TEND
  26 C
  27       DATA COEF/8.3333333333333334D-02,-2.7777777777777778D-03,
  28      .          7.9365079365079365D-04,-5.9523809523809524D-04,
  29      .          8.4175084175084175D-04,-1.9175269175269175D-03,
  30      .          6.4102564102564103D-03,-2.9550653594771242D-02,
  31      .          0.1796443723688306D0  ,-0.6962161084529506D0  /
  32       DATA PI/  3.141592653589793D0/
  33       DATA ZEROD/0.D0/,HALFD/0.5D0/,ONED/1.D0/,TWOD/2.D0/,TEND/10.D0/
  34 C
  35       X=Y
  36       AFSPL=ONED
  37       N=INT(TEND-Y)
  38       DO 10 I=0,N
  39         AFSPL=AFSPL*X
  40         X=X+ONED
  41 10    CONTINUE
  42       R=(X-HALFD)* LOG(X)-X+HALFD* LOG(TWOD*PI)
  43       S=X
  44       T=ZEROD
  45       DO 20 I=1,10
  46         T=T+COEF(I)/S
  47         S=S*X**2
  48 20    CONTINUE
  49       GAMFUN = EXP(R+T)/AFSPL
  50       END