Doxygen Source Code Documentation
cdf_01.c File Reference
#include "cdflib.h"Go to the source code of this file.
Defines | |
| #define | hln2pi 0.91893853320467274178e0 |
Functions | |
| double | alngam (double *x) |
Define Documentation
|
|
|
Function Documentation
|
|
Definition at line 2 of file cdf_01.c. References alngam(), devlpl(), fifidint(), i, and offset. Referenced by alngam(), cumchn(), and cumfnc().
00036 {
00037 #define hln2pi 0.91893853320467274178e0
00038 static double coef[5] = {
00039 0.83333333333333023564e-1,-0.27777777768818808e-2,0.79365006754279e-3,
00040 -0.594997310889e-3,0.8065880899e-3
00041 };
00042 static double scoefd[4] = {
00043 0.62003838007126989331e2,0.9822521104713994894e1,-0.8906016659497461257e1,
00044 0.1000000000000000000e1
00045 };
00046 static double scoefn[9] = {
00047 0.62003838007127258804e2,0.36036772530024836321e2,0.20782472531792126786e2,
00048 0.6338067999387272343e1,0.215994312846059073e1,0.3980671310203570498e0,
00049 0.1093115956710439502e0,0.92381945590275995e-2,0.29737866448101651e-2
00050 };
00051 static int K1 = 9;
00052 static int K3 = 4;
00053 static int K5 = 5;
00054 static double alngam,offset,prod,xx;
00055 static int i,n;
00056 static double T2,T4,T6;
00057 /*
00058 ..
00059 .. Executable Statements ..
00060 */
00061 if(!(*x <= 6.0e0)) goto S70;
00062 prod = 1.0e0;
00063 xx = *x;
00064 if(!(*x > 3.0e0)) goto S30;
00065 S10:
00066 if(!(xx > 3.0e0)) goto S20;
00067 xx -= 1.0e0;
00068 prod *= xx;
00069 goto S10;
00070 S30:
00071 S20:
00072 if(!(*x < 2.0e0)) goto S60;
00073 S40:
00074 if(!(xx < 2.0e0)) goto S50;
00075 prod /= xx;
00076 xx += 1.0e0;
00077 goto S40;
00078 S60:
00079 S50:
00080 T2 = xx-2.0e0;
00081 T4 = xx-2.0e0;
00082 alngam = devlpl(scoefn,&K1,&T2)/devlpl(scoefd,&K3,&T4);
00083 /*
00084 COMPUTE RATIONAL APPROXIMATION TO GAMMA(X)
00085 */
00086 alngam *= prod;
00087 alngam = log(alngam);
00088 goto S110;
00089 S70:
00090 offset = hln2pi;
00091 /*
00092 IF NECESSARY MAKE X AT LEAST 12 AND CARRY CORRECTION IN OFFSET
00093 */
00094 n = fifidint(12.0e0-*x);
00095 if(!(n > 0)) goto S90;
00096 prod = 1.0e0;
00097 for(i=1; i<=n; i++) prod *= (*x+(double)(i-1));
00098 offset -= log(prod);
00099 xx = *x+(double)n;
00100 goto S100;
00101 S90:
00102 xx = *x;
00103 S100:
00104 /*
00105 COMPUTE POWER SERIES
00106 */
00107 T6 = 1.0e0/pow(xx,2.0);
00108 alngam = devlpl(coef,&K5,&T6)/xx;
00109 alngam += (offset+(xx-0.5e0)*log(xx)-xx);
00110 S110:
00111 return alngam;
00112 #undef hln2pi
00113 } /* END */
|
