forester/archive/RIO/others/hmmer/testsuite/fitting_test.c

   1 /* fitting_test.c
   2  * 17 June 1997 (see notebook)
   3  */
   4
   5 #include <stdio.h>
   6 #include <stdlib.h>
   7 #include <string.h>
   8 #include <math.h>
   9 #include <time.h>
  10 #include <assert.h>
  11
  12 #include "structs.h"
  13 #include "funcs.h"
  14 #include "squid.h"
  15
  16 #include "globals.h"
  17
  18 int
  19 main(int argc, char **argv)
  20 {
  21   int   n;                      /* number of EVD samples */
  22   float p1, p2;
  23   struct histogram_s *histog;
  24   int   i,j;
  25   float x;
  26   int   seed;
  27   int   do_evd, set, fit_evd, show_hist;
  28
  29   p1 = atof(argv[1]);           /* mu or mean */
  30   p2 = atof(argv[2]);           /* lambda or sd */
  31   n  = atoi(argv[3]);           /* # of histograms */
  32   do_evd = atoi(argv[4]);       /* 1 to sample EVD; 0 to sample Gaussian */
  33   set = atoi(argv[5]);          /* 1 to set instead of fit the dist  */
  34   fit_evd = atoi(argv[6]);      /* 1 to fit EVD; 0 to fit Gaussian  */
  35   show_hist = atoi(argv[7]);    /* 1 to show histogram */
  36
  37   seed = (int) time ((time_t *) NULL);
  38   sre_srandom(seed);
  39
  40   for (j = 0; j < n; j++)
  41     {
  42       histog = AllocHistogram(-200, 200, 100);
  43       for (i = 0; i < 2500; i++)
  44         {
  45           if (do_evd) x = EVDrandom(p1, p2);
  46           else        x = Gaussrandom(p1, p2);
  47
  48           assert(x > -100.);
  49           assert(x < 100.);
  50           AddToHistogram(histog, x);
  51         }
  52
  53       if (set && fit_evd)
  54         ExtremeValueSetHistogram(histog, p1, p2);
  55       else if (set && !fit_evd)
  56         GaussianSetHistogram(histog, p1, p2);
  57       else if (!set && fit_evd)
  58         ExtremeValueFitHistogram(histog, 9999.);
  59       else
  60         GaussianFitHistogram(histog, 9999.);
  61
  62       printf("%f\n", histog->chip);
  63
  64       if (show_hist)
  65         PrintASCIIHistogram(stdout, histog);
  66
  67       FreeHistogram(histog);
  68     }
  69
  70   return 0;
  71 }