Doxygen Source Code Documentation
cdf_13.c File Reference
#include "cdflib.h"Go to the source code of this file.
Functions | |
| double | bup (double *a, double *b, double *x, double *y, int *n, double *eps) |
Function Documentation
|
||||||||||||||||||||||||||||
|
Definition at line 2 of file cdf_13.c. References a, brcmp1(), bup(), exparg(), i, l, and r. Referenced by bratio(), and bup().
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 */
|
