|
Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
cdf_19.c File Reference#include "cdflib.h"
Go to the source code of this file.
|
Defines |
#define | tent4 1.0e4 |
#define | tol (1.0e-8) |
#define | atol (1.0e-50) |
#define | zero (1.0e-300) |
#define | one (1.0e0-1.0e-16) |
#define | inf 1.0e300 |
Functions |
void | cdffnc (int *which, double *p, double *q, double *f, double *dfn, double *dfd, double *phonc, int *status, double *bound) |
Define Documentation
#define one (1.0e0-1.0e-16)
|
|
Function Documentation
void cdffnc |
( |
int * |
which, |
|
|
double * |
p, |
|
|
double * |
q, |
|
|
double * |
f, |
|
|
double * |
dfn, |
|
|
double * |
dfd, |
|
|
double * |
phonc, |
|
|
int * |
status, |
|
|
double * |
bound |
|
) |
|
|
|
Definition at line 2 of file cdf_19.c.
References cumfnc(), dinvr(), dstinv(), p, and q.
00025 : 1..5
00026 iwhich = 1 : Calculate P and Q from F,DFN,DFD and PNONC
00027 iwhich = 2 : Calculate F from P,Q,DFN,DFD and PNONC
00028 iwhich = 3 : Calculate DFN from P,Q,F,DFD and PNONC
00029 iwhich = 4 : Calculate DFD from P,Q,F,DFN and PNONC
00030 iwhich = 5 : Calculate PNONC from P,Q,F,DFN and DFD
00031
00032 P <--> The integral from 0 to F of the non-central f-density.
00033 Input range: [0,1-1E-16).
00034
00035 Q <--> 1-P.
00036 Q is not used by this subroutine and is only included
00037 for similarity with other cdf* routines.
00038
00039 F <--> Upper limit of integration of the non-central f-density.
00040 Input range: [0, +infinity).
00041 Search range: [0,1E300]
00042
00043 DFN < --> Degrees of freedom of the numerator sum of squares.
00044 Input range: (0, +infinity).
00045 Search range: [ 1E-300, 1E300]
00046
00047 DFD < --> Degrees of freedom of the denominator sum of squares.
00048 Must be in range: (0, +infinity).
00049 Input range: (0, +infinity).
00050 Search range: [ 1E-300, 1E300]
00051
00052 PNONC <-> The non-centrality parameter
00053 Input range: [0,infinity)
00054 Search range: [0,1E4]
00055
00056 STATUS <-- 0 if calculation completed correctly
00057 -I if input parameter number I is out of range
00058 1 if answer appears to be lower than lowest
00059 search bound
00060 2 if answer appears to be higher than greatest
00061 search bound
00062 3 if P + Q .ne. 1
00063
00064 BOUND <-- Undefined if STATUS is 0
00065
00066 Bound exceeded by parameter number I if STATUS
00067 is negative.
00068
00069 Lower search bound if STATUS is 1.
00070
00071 Upper search bound if STATUS is 2.
00072
00073
00074 Method
00075
00076
00077 Formula 26.6.20 of Abramowitz and Stegun, Handbook of
00078 Mathematical Functions (1966) is used to compute the cumulative
00079 distribution function.
00080
00081 Computation of other parameters involve a seach for a value that
00082 produces the desired value of P. The search relies on the
00083 monotinicity of P with the other parameter.
00084
00085 WARNING
00086
00087 The computation time required for this routine is proportional
00088 to the noncentrality parameter (PNONC). Very large values of
00089 this parameter can consume immense computer resources. This is
00090 why the search range is bounded by 10,000.
00091
00092 WARNING
00093
00094 The value of the cumulative noncentral F distribution is not
00095 necessarily monotone in either degrees of freedom. There thus
00096 may be two values that provide a given CDF value. This routine
00097 assumes monotonicity and will find an arbitrary one of the two
00098 values.
00099
00100 **********************************************************************/
00101 {
00102 #define tent4 1.0e4
00103 #define tol (1.0e-8)
00104 #define atol (1.0e-50)
00105 #define zero (1.0e-300)
00106 #define one (1.0e0-1.0e-16)
00107 #define inf 1.0e300
00108 static double K1 = 0.0e0;
00109 static double K3 = 0.5e0;
00110 static double K4 = 5.0e0;
00111 static double fx,cum,ccum;
00112 static unsigned long qhi,qleft;
00113 static double T2,T5,T6,T7,T8,T9,T10,T11,T12,T13,T14,T15,T16,T17;
00114
00115
00116
00117
00118
00119
00120
00121 if(!(*which < 1 || *which > 5)) goto S30;
00122 if(!(*which < 1)) goto S10;
00123 *bound = 1.0e0;
00124 goto S20;
00125 S10:
00126 *bound = 5.0e0;
00127 S20:
00128 *status = -1;
00129 return;
00130 S30:
00131 if(*which == 1) goto S70;
00132
00133
00134
00135 if(!(*p < 0.0e0 || *p > one)) goto S60;
00136 if(!(*p < 0.0e0)) goto S40;
00137 *bound = 0.0e0;
00138 goto S50;
00139 S40:
00140 *bound = one;
00141 S50:
00142 *status = -2;
00143 return;
00144 S70:
00145 S60:
00146 if(*which == 2) goto S90;
00147
00148
00149
00150 if(!(*f < 0.0e0)) goto S80;
00151 *bound = 0.0e0;
00152 *status = -4;
00153 return;
00154 S90:
00155 S80:
00156 if(*which == 3) goto S110;
00157
00158
00159
00160 if(!(*dfn <= 0.0e0)) goto S100;
00161 *bound = 0.0e0;
00162 *status = -5;
00163 return;
00164 S110:
00165 S100:
00166 if(*which == 4) goto S130;
00167
00168
00169
00170 if(!(*dfd <= 0.0e0)) goto S120;
00171 *bound = 0.0e0;
00172 *status = -6;
00173 return;
00174 S130:
00175 S120:
00176 if(*which == 5) goto S150;
00177
00178
00179
00180 if(!(*phonc < 0.0e0)) goto S140;
00181 *bound = 0.0e0;
00182 *status = -7;
00183 return;
00184 S150:
00185 S140:
00186
00187
00188
00189 if(1 == *which) {
00190
00191
00192
00193 cumfnc(f,dfn,dfd,phonc,p,q);
00194 *status = 0;
00195 }
00196 else if(2 == *which) {
00197
00198
00199
00200 *f = 5.0e0;
00201 T2 = inf;
00202 T5 = atol;
00203 T6 = tol;
00204 dstinv(&K1,&T2,&K3,&K3,&K4,&T5,&T6);
00205 *status = 0;
00206 dinvr(status,f,&fx,&qleft,&qhi);
00207 S160:
00208 if(!(*status == 1)) goto S170;
00209 cumfnc(f,dfn,dfd,phonc,&cum,&ccum);
00210 fx = cum-*p;
00211 dinvr(status,f,&fx,&qleft,&qhi);
00212 goto S160;
00213 S170:
00214 if(!(*status == -1)) goto S200;
00215 if(!qleft) goto S180;
00216 *status = 1;
00217 *bound = 0.0e0;
00218 goto S190;
00219 S180:
00220 *status = 2;
00221 *bound = inf;
00222 S200:
00223 S190:
00224 ;
00225 }
00226 else if(3 == *which) {
00227
00228
00229
00230 *dfn = 5.0e0;
00231 T7 = zero;
00232 T8 = inf;
00233 T9 = atol;
00234 T10 = tol;
00235 dstinv(&T7,&T8,&K3,&K3,&K4,&T9,&T10);
00236 *status = 0;
00237 dinvr(status,dfn,&fx,&qleft,&qhi);
00238 S210:
00239 if(!(*status == 1)) goto S220;
00240 cumfnc(f,dfn,dfd,phonc,&cum,&ccum);
00241 fx = cum-*p;
00242 dinvr(status,dfn,&fx,&qleft,&qhi);
00243 goto S210;
00244 S220:
00245 if(!(*status == -1)) goto S250;
00246 if(!qleft) goto S230;
00247 *status = 1;
00248 *bound = zero;
00249 goto S240;
00250 S230:
00251 *status = 2;
00252 *bound = inf;
00253 S250:
00254 S240:
00255 ;
00256 }
00257 else if(4 == *which) {
00258
00259
00260
00261 *dfd = 5.0e0;
00262 T11 = zero;
00263 T12 = inf;
00264 T13 = atol;
00265 T14 = tol;
00266 dstinv(&T11,&T12,&K3,&K3,&K4,&T13,&T14);
00267 *status = 0;
00268 dinvr(status,dfd,&fx,&qleft,&qhi);
00269 S260:
00270 if(!(*status == 1)) goto S270;
00271 cumfnc(f,dfn,dfd,phonc,&cum,&ccum);
00272 fx = cum-*p;
00273 dinvr(status,dfd,&fx,&qleft,&qhi);
00274 goto S260;
00275 S270:
00276 if(!(*status == -1)) goto S300;
00277 if(!qleft) goto S280;
00278 *status = 1;
00279 *bound = zero;
00280 goto S290;
00281 S280:
00282 *status = 2;
00283 *bound = inf;
00284 S300:
00285 S290:
00286 ;
00287 }
00288 else if(5 == *which) {
00289
00290
00291
00292 *phonc = 5.0e0;
00293 T15 = tent4;
00294 T16 = atol;
00295 T17 = tol;
00296 dstinv(&K1,&T15,&K3,&K3,&K4,&T16,&T17);
00297 *status = 0;
00298 dinvr(status,phonc,&fx,&qleft,&qhi);
00299 S310:
00300 if(!(*status == 1)) goto S320;
00301 cumfnc(f,dfn,dfd,phonc,&cum,&ccum);
00302 fx = cum-*p;
00303 dinvr(status,phonc,&fx,&qleft,&qhi);
00304 goto S310;
00305 S320:
00306 if(!(*status == -1)) goto S350;
00307 if(!qleft) goto S330;
00308 *status = 1;
00309 *bound = 0.0e0;
00310 goto S340;
00311 S330:
00312 *status = 2;
00313 *bound = tent4;
00314 S340:
00315 ;
00316 }
00317 S350:
00318 return;
00319 #undef tent4
00320 #undef tol
00321 #undef atol
00322 #undef zero
00323 #undef one
00324 #undef inf
00325 }
|
|