|
|
Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
cdf_24.c File Reference#include "cdflib.h"
Go to the source code of this file.
|
Defines |
| #define | tol (1.0e-8) |
| #define | atol (1.0e-50) |
| #define | zero (1.0e-300) |
| #define | inf 1.0e300 |
| #define | maxdf 1.0e10 |
Functions |
| void | cdft (int *which, double *p, double *q, double *t, double *df, int *status, double *bound) |
Define Documentation
Function Documentation
| void cdft |
( |
int * |
which, |
|
|
double * |
p, |
|
|
double * |
q, |
|
|
double * |
t, |
|
|
double * |
df, |
|
|
int * |
status, |
|
|
double * |
bound |
|
) |
|
|
|
|
Definition at line 2 of file cdf_24.c.
References cumt(), dinvr(), dstinv(), dt1(), p, q, and spmpar().
00025 : 1..3
00026 iwhich = 1 : Calculate P and Q from T and DF
00027 iwhich = 2 : Calculate T from P,Q and DF
00028 iwhich = 3 : Calculate DF from P,Q and T
00029
00030 P <--> The integral from -infinity to t of the t-density.
00031 Input range: (0,1].
00032
00033 Q <--> 1-P.
00034 Input range: (0, 1].
00035 P + Q = 1.0.
00036
00037 T <--> Upper limit of integration of the t-density.
00038 Input range: ( -infinity, +infinity).
00039 Search range: [ -1E300, 1E300 ]
00040
00041 DF <--> Degrees of freedom of the t-distribution.
00042 Input range: (0 , +infinity).
00043 Search range: [1e-300, 1E10]
00044
00045 STATUS <-- 0 if calculation completed correctly
00046 -I if input parameter number I is out of range
00047 1 if answer appears to be lower than lowest
00048 search bound
00049 2 if answer appears to be higher than greatest
00050 search bound
00051 3 if P + Q .ne. 1
00052
00053 BOUND <-- Undefined if STATUS is 0
00054
00055 Bound exceeded by parameter number I if STATUS
00056 is negative.
00057
00058 Lower search bound if STATUS is 1.
00059
00060 Upper search bound if STATUS is 2.
00061
00062
00063 Method
00064
00065
00066 Formula 26.5.27 of Abramowitz and Stegun, Handbook of
00067 Mathematical Functions (1966) is used to reduce the computation
00068 of the cumulative distribution function to that of an incomplete
00069 beta.
00070
00071 Computation of other parameters involve a seach for a value that
00072 produces the desired value of P. The search relies on the
00073 monotinicity of P with the other parameter.
00074
00075 **********************************************************************/
00076 {
00077 #define tol (1.0e-8)
00078 #define atol (1.0e-50)
00079 #define zero (1.0e-300)
00080 #define inf 1.0e300
00081 #define maxdf 1.0e10
00082 static int K1 = 1;
00083 static double K4 = 0.5e0;
00084 static double K5 = 5.0e0;
00085 static double fx,cum,ccum,pq;
00086 static unsigned long qhi,qleft,qporq;
00087 static double T2,T3,T6,T7,T8,T9,T10,T11;
00088
00089
00090
00091
00092
00093
00094
00095 if(!(*which < 1 || *which > 3)) goto S30;
00096 if(!(*which < 1)) goto S10;
00097 *bound = 1.0e0;
00098 goto S20;
00099 S10:
00100 *bound = 3.0e0;
00101 S20:
00102 *status = -1;
00103 return;
00104 S30:
00105 if(*which == 1) goto S70;
00106
00107
00108
00109 if(!(*p <= 0.0e0 || *p > 1.0e0)) goto S60;
00110 if(!(*p <= 0.0e0)) goto S40;
00111 *bound = 0.0e0;
00112 goto S50;
00113 S40:
00114 *bound = 1.0e0;
00115 S50:
00116 *status = -2;
00117 return;
00118 S70:
00119 S60:
00120 if(*which == 1) goto S110;
00121
00122
00123
00124 if(!(*q <= 0.0e0 || *q > 1.0e0)) goto S100;
00125 if(!(*q <= 0.0e0)) goto S80;
00126 *bound = 0.0e0;
00127 goto S90;
00128 S80:
00129 *bound = 1.0e0;
00130 S90:
00131 *status = -3;
00132 return;
00133 S110:
00134 S100:
00135 if(*which == 3) goto S130;
00136
00137
00138
00139 if(!(*df <= 0.0e0)) goto S120;
00140 *bound = 0.0e0;
00141 *status = -5;
00142 return;
00143 S130:
00144 S120:
00145 if(*which == 1) goto S170;
00146
00147
00148
00149 pq = *p+*q;
00150 if(!(fabs(pq-0.5e0-0.5e0) > 3.0e0*spmpar(&K1))) goto S160;
00151 if(!(pq < 0.0e0)) goto S140;
00152 *bound = 0.0e0;
00153 goto S150;
00154 S140:
00155 *bound = 1.0e0;
00156 S150:
00157 *status = 3;
00158 return;
00159 S170:
00160 S160:
00161 if(!(*which == 1)) qporq = *p <= *q;
00162
00163
00164
00165
00166 if(1 == *which) {
00167
00168
00169
00170 cumt(t,df,p,q);
00171 *status = 0;
00172 }
00173 else if(2 == *which) {
00174
00175
00176
00177
00178 *t = dt1(p,q,df);
00179 T2 = -inf;
00180 T3 = inf;
00181 T6 = atol;
00182 T7 = tol;
00183 dstinv(&T2,&T3,&K4,&K4,&K5,&T6,&T7);
00184 *status = 0;
00185 dinvr(status,t,&fx,&qleft,&qhi);
00186 S180:
00187 if(!(*status == 1)) goto S210;
00188 cumt(t,df,&cum,&ccum);
00189 if(!qporq) goto S190;
00190 fx = cum-*p;
00191 goto S200;
00192 S190:
00193 fx = ccum-*q;
00194 S200:
00195 dinvr(status,t,&fx,&qleft,&qhi);
00196 goto S180;
00197 S210:
00198 if(!(*status == -1)) goto S240;
00199 if(!qleft) goto S220;
00200 *status = 1;
00201 *bound = -inf;
00202 goto S230;
00203 S220:
00204 *status = 2;
00205 *bound = inf;
00206 S240:
00207 S230:
00208 ;
00209 }
00210 else if(3 == *which) {
00211
00212
00213
00214 *df = 5.0e0;
00215 T8 = zero;
00216 T9 = maxdf;
00217 T10 = atol;
00218 T11 = tol;
00219 dstinv(&T8,&T9,&K4,&K4,&K5,&T10,&T11);
00220 *status = 0;
00221 dinvr(status,df,&fx,&qleft,&qhi);
00222 S250:
00223 if(!(*status == 1)) goto S280;
00224 cumt(t,df,&cum,&ccum);
00225 if(!qporq) goto S260;
00226 fx = cum-*p;
00227 goto S270;
00228 S260:
00229 fx = ccum-*q;
00230 S270:
00231 dinvr(status,df,&fx,&qleft,&qhi);
00232 goto S250;
00233 S280:
00234 if(!(*status == -1)) goto S310;
00235 if(!qleft) goto S290;
00236 *status = 1;
00237 *bound = zero;
00238 goto S300;
00239 S290:
00240 *status = 2;
00241 *bound = maxdf;
00242 S300:
00243 ;
00244 }
00245 S310:
00246 return;
00247 #undef tol
00248 #undef atol
00249 #undef zero
00250 #undef inf
00251 #undef maxdf
00252 }
|
|
