C      ALGORITHM 344, COLLECTED ALGORITHMS FROM ACM.
C      THIS WORK PUBLISHED IN COMMUNICATIONS OF THE ACM
C      VOL. 12, NO. 1, January, 1969, PP.37--38.
#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create the files:
#	Fortran/
#	Fortran/Sp/
#	Fortran/Sp/Drivers/
#	Fortran/Sp/Drivers/Makefile
#	Fortran/Sp/Drivers/driver.f
#	Fortran/Sp/Drivers/res
#	Fortran/Sp/Src/
#	Fortran/Sp/Src/src.f
# This archive created: Wed Jan 18 20:30:20 2006
export PATH; PATH=/bin:$PATH
if test ! -d 'Fortran'
then
	mkdir 'Fortran'
fi
cd 'Fortran'
if test ! -d 'Sp'
then
	mkdir 'Sp'
fi
cd 'Sp'
if test ! -d 'Drivers'
then
	mkdir 'Drivers'
fi
cd 'Drivers'
if test -f 'Makefile'
then
	echo shar: will not over-write existing file "'Makefile'"
else
cat << "SHAR_EOF" > 'Makefile'
all: Res

src.o: src.f
	$(F77) $(F77OPTS) -c src.f

driver.o: driver.f
	$(F77) $(F77OPTS) -c driver.f

DRIVERS= driver
RESULTS= Res

Objs1= driver.o src.o
driver: $(Objs1)
	$(F77) $(F77OPTS) -o driver $(Objs1) $(SRCLIBS)
Res: driver 
	./driver >Res

diffres:Res res
	echo "Differences in results from driver"
	$(DIFF) Res res

clean: 
	rm -rf *.o $(DRIVERS) $(CLEANUP) $(RESULTS)
SHAR_EOF
fi # end of overwriting check
if test -f 'driver.f'
then
	echo shar: will not over-write existing file "'driver.f'"
else
cat << "SHAR_EOF" > 'driver.f'
      program main

c***********************************************************************
c
cc TOMS344_PRB tests TTEST
c
      implicit none

      write ( *, '(a)' ) ' '
      write ( *, '(a)' ) 'TOMS344_PRB:'
      write ( *, '(a)' ) '  Test ACM algorithm 344, which evaluates'
      write ( *, '(a)' ) '  Student''s T distribution.'
 
      call test01

      write ( *, '(a)' ) ' '
      write ( *, '(a)' ) 'TOMS344_PRB'
      write ( *, '(a)' ) '  Normal end of execution.'

      stop
      end
      subroutine test01

c*******************************************************************************
c
cc TEST01 compares TTEST versus values from STUDENT_CDF_VALUES.
c
c  Modified:
c
c    06 January 2006
c
      implicit none

      real c
      integer df
      real fx
      real fx2
      integer kerr
      integer n_data
      real x

      write ( *, '(a)' ) ' '
      write ( *, '(a)' ) 'TEST140:'
      write ( *, '(a)' ) '  STUDENT_CDF_VALUES returns values of '
      write ( *, '(a)' ) '  the Student T CDF.'
      write ( *, '(a)' ) ' '
      write ( *, '(a)' ) 
     &  '      DF        X         CDF                       CDF'
      write ( *, '(a)' ) 
     &  '                          exact                     computed'
      write ( *, '(a)' ) ' '

      n_data = 0

 10   continue

      call student_cdf_values ( n_data, c, x, fx )

      if ( n_data == 0 ) then
        go to 20
      end if

      df = int ( c )

      if ( c .ne. real ( df ) ) then
        go to 10
      end if

      call ttest ( x, df, fx2, kerr )
c
c  Transform FX.  TTEST computes the two-tailed distribution, and uses the
c  complementary form of the integral.
c
      fx = 2.0E+00 * ( 1.0E+00 - fx )

      write ( *, '(2x,i6,2x,f10.4,2x,g24.16,2x,g24.16)' ) df, x, fx, fx2

      go to 10

20    continue

      return
      end
      subroutine student_cdf_values ( n_data, c, x, fx )

c*******************************************************************************
c
cc STUDENT_CDF_VALUES returns some values of the Student CDF.
c
c  Discussion:
c
c    In Mathematica, the function can be evaluated by:
c
c      Needs["Statistics`ContinuousDistributions`"]
c      dist = StudentTDistribution [ c ]
c      CDF [ dist, x ]
c
c  Modified:
c
c    06 January 2005
c
c  Author:
c
c    John Burkardt
c
c  Reference:
c
c    Milton Abramowitz and Irene Stegun,
c    Handbook of Mathematical Functions,
c    US Department of Commerce, 1964.
c
c    Stephen Wolfram,
c    The Mathematica Book,
c    Fourth Edition,
c    Wolfram Media / Cambridge University Press, 1999.
c
c  Parameters:
c
c    Input/output, integer N_DATA.  The user sets N_DATA to 0 before the
c    first call.  On each call, the routine increments N_DATA by 1, and
c    returns the corresponding data; when there is no more data, the
c    output value of N_DATA will be 0 again.
c
c    Output, real C, is usually called the number of
c    degrees of freedom of the distribution.  C is typically an
c    integer, but that is not essential.  It is required that
c    C be strictly positive.
c
c    Output, real X, the argument of the function.
c
c    Output, real FX, the value of the function.
c
      implicit none

      integer n_max
      parameter ( n_max = 13 )

      real c
      real c_vec(n_max)
      real fx
      real fx_vec(n_max)
      integer n_data
      real x
      real x_vec(n_max)

      save c_vec
      save fx_vec
      save x_vec

      data c_vec /
     &  1.0E+00,
     &  2.0E+00,
     &  3.0E+00,
     &  4.0E+00,
     &  5.0E+00,
     &  2.0E+00,
     &  5.0E+00,
     &  2.0E+00,
     &  5.0E+00,
     &  2.0E+00,
     &  3.0E+00,
     &  4.0E+00,
     &  5.0E+00 /
      data fx_vec /
     &  0.6000231200328521E+00,
     &  0.6001080279134390E+00,
     &  0.6001150934648930E+00,
     &  0.6000995134721354E+00,
     &  0.5999341989834830E+00,
     &  0.7498859393137811E+00,
     &  0.7500879487671045E+00,
     &  0.9500004222186464E+00,
     &  0.9499969138365968E+00,
     &  0.9900012348724744E+00,
     &  0.9900017619355059E+00,
     &  0.9900004567580596E+00,
     &  0.9900007637471291E+00 /
      data x_vec /
     &  0.325E+00,
     &  0.289E+00,
     &  0.277E+00,
     &  0.271E+00,
     &  0.267E+00,
     &  0.816E+00,
     &  0.727E+00,
     &  2.920E+00,
     &  2.015E+00,
     &  6.965E+00,
     &  4.541E+00,
     &  3.747E+00,
     &  3.365E+00 /

      if ( n_data < 0 ) then
        n_data = 0
      end if

      n_data = n_data + 1

      if ( n_max < n_data ) then
        n_data = 0
        c = 0.0E+00
        x = 0.0E+00
        fx = 0.0E+00
      else
        c = c_vec(n_data)
        x = x_vec(n_data)
        fx = fx_vec(n_data)
      end if

      return
      end
SHAR_EOF
fi # end of overwriting check
if test -f 'res'
then
	echo shar: will not over-write existing file "'res'"
else
cat << "SHAR_EOF" > 'res'
 
TOMS344_PRB:
  Test ACM algorithm 344, which evaluates
  Student's T distribution.
 
TEST140:
  STUDENT_CDF_VALUES returns values of 
  the Student T CDF.
 
      DF        X         CDF                       CDF
                          exact                     computed
 
       1      0.3250    0.7999538183212280        0.7999537587165833    
       2      0.2890    0.7997839450836182        0.7997839450836182    
       3      0.2770    0.7997697591781616        0.7997698187828064    
       4      0.2710    0.7998009920120239        0.7998009920120239    
       5      0.2670    0.8001315593719482        0.8001316189765930    
       2      0.8160    0.5002281665802002        0.5002281069755554    
       5      0.7270    0.4998240470886230        0.4998240470886230    
       2      2.9200    0.9999918937683105E-01    0.9999912977218628E-01
       5      2.0150    0.1000062227249146        0.1000061407685280    
       2      6.9650    0.1999747753143311E-01    0.1999753713607788E-01
       3      4.5410    0.1999652385711670E-01    0.1999655365943909E-01
       4      3.7470    0.1999902725219727E-01    0.1999907195568085E-01
       5      3.3650    0.1999843120574951E-01    0.1999859884381294E-01
 
TOMS344_PRB
  Normal end of execution.
SHAR_EOF
fi # end of overwriting check
cd ..
if test ! -d 'Src'
then
	mkdir 'Src'
fi
cd 'Src'
if test -f 'src.f'
then
	echo shar: will not over-write existing file "'src.f'"
else
cat << "SHAR_EOF" > 'src.f'
      SUBROUTINE TTEST ( T, DF, ANS, KERR )
C
      IMPLICIT NONE

      REAL ANS
      REAL D1
      REAL D2
      INTEGER DF
      REAL F1
      REAL F2
      INTEGER I
      INTEGER KERR
      INTEGER N
      REAL T
      REAL T1
      REAL T2

      SAVE D1

      DATA D1 / 0.63661977236758134E+00 /
C
C  0.63661977236758134... = 2 / PI
C
      KERR = 0
C
      IF ( DF .GT. 0 ) GO TO 1
C
C  ERROR RETURN IF DF NOT POSITIVE.
C
      KERR = 1
      ANS = 0.0E+00
      RETURN
C
C  BEGIN COMPUTATION OF SERIES.
C
1     T2 = T * T / FLOAT ( DF )
      T1 = SQRT ( T2 )
      T2 = 1.0E+00 / ( 1.0E+00 + T2 )
C
      IF ( ( DF / 2 ) * 2 .EQ. DF ) GO TO 5
C
C  DF IS AN ODD INTEGER.
C
      ANS = 1.0E+00 - D1 * ATAN ( T1 )
C
      IF ( DF .EQ. 1 ) GO TO 4
C
      D2 = D1 * T1 * T2
      ANS = ANS - D2
C
      IF ( DF .EQ. 3 ) GO TO 4
C
      F1 = 0.0E+00
2     N = ( DF - 2 ) / 2
      DO 3 I = 1, N
        F2 = 2.0E+00 * FLOAT ( I ) - F1
        D2 = D2 * T2 * F2 / ( F2 + 1.0E+00 )
        ANS = ANS - D2
3     CONTINUE
C
C  COMMON RETURN AFTER COMPUTATION.
C
4      IF ( ANS .LT. 0.0E+00 ) ANS = 0.0E+00
       RETURN
C
C  DF IS AN EVEN INTEGER.
C
5      D2 = T1 * SQRT ( T2 )
       ANS = 1.0E+00 - D2
C
       IF ( DF .EQ. 2 ) GO TO 4
C
       F1 = 1.0E+00
       GO TO 2
       END
SHAR_EOF
fi # end of overwriting check
cd ..
cd ..
cd ..
#       End of shell archive
exit 0