Time sequence and spectral analysis - Earth Dynamics: Deformations and Oscillations of the Rotating Earth - page 193

C Set initial values of even and odd sequences a and b.
AAE=1.D0-A
AAO=1.D0
BBE=1.D0
BBO=X
C Evaluate continued fraction, first by odd sequences, then by even.
DO 10 N=1,ITER
C Calculate odd a and b.
AO=BBO*AE+AAO*AO
BO=BBO*BE+AAO*BO
CF=AO/BO
C Test for convergence.
IF(CF.EQ.CFOLD) GO TO 11
C Store old CF.
CFOLD=CF
C Update odd sequence a.
IF(N.GT.1)AAO=AAO+1.D0
C Calculate even a and b.
AE=BBE*AO+AAE*AE
BE=BBE*BO+AAE*BE
CF=AE/BE
C Update even sequence a.
AAE=AAE+1.D0
C Test for convergence.
IF(CF.EQ.CFOLD) GO TO 11
C Store old CF.
CFOLD=CF

10 CONTINUE

PAUSE 'A too large or X too small.'

11 CFGI=1.D0-DEXP(-X+A*DLOG(X)-ALG)*CF

RETURN
END

The subroutine GAMI calls the subroutine SEGAMI or the subroutine CFGAMI
to compute the incomplete gamma function for all arguments for which A is posit-
ive and X is non-negative.

SUBROUTINE GAMI(GI,A,X)

C
C Subroutine for the calculation of the incomplete gamma function.
C

IMPLICIT DOUBLE PRECISION(A-H,O-Z)
C Check if X is negative or A is not positive.
IF(X.LT.0.D0.OR.A.LE.0.D0) GO TO 10
C If X is less than A+1, use series calculation,
C otherwise use continued fraction expansion.
IF(X.LT.A+1.D0) CALL SEGAMI(GI,A,X)
IF(X.GE.A+1.D0) CALL CFGAMI(GI,A,X)
GO TO 11

10 PAUSE 'X is negative or A is not positive.'

11 CONTINUE

RETURN
END

Evaluation of the cumulative distribution function F ( z ) given by expression
(2.397) requires the equivalent degrees of freedom of each individual spectrum
entering the product. For 75% overlap, the equivalent number of degrees of