home *** CD-ROM | disk | FTP | other *** search
/ QBasic & Borland Pascal & C / Delphi5.iso / C / Samples / C-SSP.ARJ / STHYGD.C < prev    next >
Encoding:
Text File  |  1984-08-30  |  1.7 KB  |  81 lines
  1.    sthygd(k,n,id,iix,df,cf)
  2.  
  3.       /* this function computes the probability function and cumulative */
  4.       /* function for hypergeometric distribution of the form           */
  5.       /* df = c(id,iix)*c(k-id,n-iix)/c(k,n).                           */
  6.  
  7.       int id,iix,k,n;
  8.       double *cf,*df;
  9.  
  10.    {
  11.  
  12.       int i,idx,ix,ix1,kd,kdnx,kn,n1,nx,nx1;
  13.       double sd,sdx,fi,cd,ckd,sp,ckdlog,cklog,ck,
  14.              skdlog,sknxlg,sklog,sknlog;
  15.       extern double exp(),log();
  16.  
  17.       ix = 0;
  18.       *cf = 0.;
  19.       while((ix-iix) <= 0)
  20.       {
  21.        idx = id - ix;
  22.        sd = 1.;
  23.        sdx = 1.;
  24.        ix1 = ix + 1;
  25.  
  26.        for(i = ix1; i <= id; i++)
  27.        {
  28.         fi = i;
  29.         sd = sd * fi;
  30.        }
  31.  
  32.        for(i = 1; i <= idx; i++)
  33.        {
  34.         fi = i;
  35.         sdx = sdx * fi;
  36.        }
  37.        cd = sd / sdx;
  38.        kd = k - id;
  39.        nx = n - ix;
  40.        nx1 = nx + 1;
  41.        skdlog = 0.;
  42.        sknxlg = 0.;
  43.        kdnx = kd - nx;
  44.  
  45.        for(i = nx1; i <= kd; i++)
  46.        {
  47.         fi = i;
  48.         skdlog = skdlog + log(fi);
  49.        }
  50.  
  51.        for(i = 1; i <= kdnx; i++)
  52.        {
  53.         fi = i;
  54.         sknxlg = sknxlg + log(fi);
  55.        }
  56.        ckdlog = skdlog - sknxlg;
  57.        ckd = exp(ckdlog);
  58.        sp = cd * ckd;
  59.        *cf = *cf + sp;
  60.        ix++;
  61.       }
  62.       sklog = 0.;
  63.       sknlog = 0.;
  64.       n1 = n + 1;
  65.       kn = k - n;
  66.       for(i = n1; i <= k; i++)
  67.       {
  68.        fi = i;
  69.        sklog = sklog + log(fi);
  70.       }
  71.       for(i = 1; i <= kn; i++)
  72.       {
  73.        fi = i;
  74.        sknlog = sknlog + log(fi);
  75.       }
  76.       cklog = sklog - sknlog;
  77.       ck = exp(cklog);
  78.       *cf = *cf / ck;
  79.       *df = sp / ck;
  80.    }
  81.