Doxygen Source Code Documentation

Main Page Alphabetical List Data Structures File List Data Fields Globals Search
cdf_13.c

00001 #include "cdflib.h"
00002 double bup(double *a,double *b,double *x,double *y,int *n,double *eps)
00003 /*
00004 -----------------------------------------------------------------------
00005      EVALUATION OF IX(A,B) - IX(A+N,B) WHERE N IS A POSITIVE INTEGER.
00006      EPS IS THE TOLERANCE USED.
00007 -----------------------------------------------------------------------
00008 */
00009 {
00010 static int K1 = 1;
00011 static int K2 = 0;
00012 static double bup,ap1,apb,d,l,r,t,w;
00013 static int i,k,kp1,mu,nm1;
00014 /*
00015      ..
00016      .. Executable Statements ..
00017 */
00018 /*
00019           OBTAIN THE SCALING FACTOR EXP(-MU) AND
00020              EXP(MU)*(X**A*Y**B/BETA(A,B))/A
00021 */
00022     apb = *a+*b;
00023     ap1 = *a+1.0e0;
00024     mu = 0;
00025     d = 1.0e0;
00026     if(*n == 1 || *a < 1.0e0) goto S10;
00027     if(apb < 1.1e0*ap1) goto S10;
00028     mu = fabs(exparg(&K1));
00029     k = exparg(&K2);
00030     if(k < mu) mu = k;
00031     t = mu;
00032     d = exp(-t);
00033 S10:
00034     bup = brcmp1(&mu,a,b,x,y)/ *a;
00035     if(*n == 1 || bup == 0.0e0) return bup;
00036     nm1 = *n-1;
00037     w = d;
00038 /*
00039           LET K BE THE INDEX OF THE MAXIMUM TERM
00040 */
00041     k = 0;
00042     if(*b <= 1.0e0) goto S50;
00043     if(*y > 1.e-4) goto S20;
00044     k = nm1;
00045     goto S30;
00046 S20:
00047     r = (*b-1.0e0)**x/ *y-*a;
00048     if(r < 1.0e0) goto S50;
00049     k = t = nm1;
00050     if(r < t) k = r;
00051 S30:
00052 /*
00053           ADD THE INCREASING TERMS OF THE SERIES
00054 */
00055     for(i=1; i<=k; i++) {
00056         l = i-1;
00057         d = (apb+l)/(ap1+l)**x*d;
00058         w += d;
00059     }
00060     if(k == nm1) goto S70;
00061 S50:
00062 /*
00063           ADD THE REMAINING TERMS OF THE SERIES
00064 */
00065     kp1 = k+1;
00066     for(i=kp1; i<=nm1; i++) {
00067         l = i-1;
00068         d = (apb+l)/(ap1+l)**x*d;
00069         w += d;
00070         if(d <= *eps*w) goto S70;
00071     }
00072 S70:
00073 /*
00074                TERMINATE THE PROCEDURE
00075 */
00076     bup *= w;
00077     return bup;
00078 } /* END */
AFNI/NIfTI Server

Sections

Personal tools

Navigation

Doxygen Source Code Documentation

cdf_13.c
