nstats.c
上传用户:lengbin
上传日期:2010-03-31
资源大小:121k
文件大小:13k
开发平台:

C/C++

  1. /*----------------------------------------------------------------------
  2.   File    : nstats.c
  3.   Contents: management of normalization statistics
  4.   Author  : Christian Borgelt
  5.   History : 12.08.2003 file created
  6.             12.08.2004 description and parse function added
  7. ----------------------------------------------------------------------*/
  8. #include <stdio.h>
  9. #include <stdlib.h>
  10. #include <string.h>
  11. #include <float.h>
  12. #include <math.h>
  13. #include <assert.h>
  14. #include "nstats.h"
  15. /*----------------------------------------------------------------------
  16.   Preprocessor Definitions
  17. ----------------------------------------------------------------------*/
  18. #define BLKSIZE   64            /* block size for parsing */
  19. /*----------------------------------------------------------------------
  20.   Functions
  21. ----------------------------------------------------------------------*/
  22. NSTATS* nst_create (int dim)
  23. {                               /* --- create numerical statistics */
  24.   NSTATS *nst;                  /* created statistics structure */
  25.   double *p;                    /* to organize the memory */
  26.   assert(dim > 0);              /* check the function argument */
  27.   nst = (NSTATS*)malloc(sizeof(NSTATS) +(6*dim -1) *sizeof(double));
  28.   if (!nst) return NULL;        /* create a statistics structure */
  29.   nst->dim  = dim;              /* and initialize the fields */
  30.   nst->reg  = 0;
  31.   nst->offs = p = nst->facs +dim;
  32.   nst->mins = p += dim;
  33.   nst->maxs = p += dim;         /* organize the vectors */
  34.   nst->sums = p += dim;
  35.   nst->sqrs = p += dim;
  36.   while (--dim >= 0) {          /* traverse the vectors */
  37.     nst->mins[dim] = DBL_MAX; nst->maxs[dim] = -DBL_MAX;
  38.     nst->sums[dim] = nst->sqrs[dim] = nst->offs[dim] = 0;
  39.     nst->facs[dim] = 1;         /* initialize the ranges of values */
  40.   }                             /* and the aggregation variables */
  41.   return nst;                   /* return created structure */
  42. }  /* nst_create() */
  43. /*--------------------------------------------------------------------*/
  44. void nst_delete (NSTATS *nst)
  45. { free(nst); }                  /* --- delete numerical statistics */
  46. /*--------------------------------------------------------------------*/
  47. void nst_reg (NSTATS *nst, const double *vec, double weight)
  48. {                               /* --- register a data vector */
  49.   int    i;                     /* loop variable */
  50.   double *min, *max;            /* to traverse the min./max. values */
  51.   double *sum, *sqr;            /* to traverse the value sums */
  52.   double *off, *fac;            /* to traverse the offsets/scales */
  53.   double t;                     /* temporary buffer */
  54.   assert(nst && vec);           /* check the function arguments */
  55.   sum = nst->sums;              /* get the vectors for the sums */
  56.   sqr = nst->sqrs;              /* and the sums of squares */
  57.   if (!vec) {                   /* if to terminate registration */
  58.     off = nst->offs;            /* get the offsets and */
  59.     fac = nst->facs;            /* the scaling factors */
  60.     if (nst->reg <= 0)          /* if no patterns are registered */
  61.       for (i = nst->dim; --i >= 0; ) { off[i] = 0; fac[i] = 1; }
  62.     else {                      /* if patterns have been registered */
  63.       for (i = nst->dim; --i >= 0; ) {      /* traverse the vectors */
  64.         off[i] = sum[i] /nst->reg;
  65.         t      = sqr[i] -off[i] *sum[i];
  66.         fac[i] = (t > 0) ? sqrt(nst->reg /t) : 1;
  67.       }                         /* estimate the parameters */
  68.     }
  69.     if (weight < 0) {           /* if to reinitialize registration */
  70.       for (i = nst->dim; --i >= 0; )
  71.         sum[i] = sqr[i] = 0;    /* reinitialize the vectors */
  72.       nst->reg = 0;             /* and the pattern counter */
  73.     } }
  74.   else {                        /* if to register a data vector */
  75.     min = nst->mins;            /* get the minimal */
  76.     max = nst->maxs;            /* and the maximal values */
  77.     for (i = nst->dim; --i >= 0; ) {
  78.       if (vec[i] < min[i]) min[i] = vec[i];
  79.       if (vec[i] > max[i]) max[i] = vec[i];
  80.       sum[i] += vec[i];         /* update the ranges of values */
  81.       sqr[i] += vec[i] *vec[i]; /* and sum the values */
  82.     }                           /* and their squares */
  83.     nst->reg += weight;         /* count the pattern */
  84.   }
  85. }  /* nst_reg() */
  86. /*--------------------------------------------------------------------*/
  87. void nst_range (NSTATS *nst, int idx, double min, double max)
  88. {                               /* --- set range of values */
  89.   int i;                        /* loop variable */
  90.   assert(nst && (idx < nst->dim));  /* check the arguments */
  91.   if (idx < 0) { i = nst->dim; idx = 0; }
  92.   else         { i = idx +1; }  /* get index range to set */
  93.   while (--i >= idx) {          /* and traverse it */
  94.     nst->mins[i] = min;         /* set the minimal */
  95.     nst->maxs[i] = max;         /* and the maximal value */
  96.   }                             /* for all dimensions in range */
  97. }  /* nst_range() */
  98. /*--------------------------------------------------------------------*/
  99. void nst_expand (NSTATS *nst, int idx, double factor)
  100. {                               /* --- expand range of values */
  101.   int    i;                     /* loop variable */
  102.   double t;                     /* change of minimal/maximal value */
  103.   assert(nst                    /* check the function arguments */
  104.      && (idx < nst->dim) && (factor >= 0));
  105.   if (idx < 0) { i = nst->dim; idx = 0; }
  106.   else         { i = idx +1; }  /* get index range to expand */
  107.   while (--i >= idx) {          /* and traverse it */
  108.     t = (nst->maxs[i] -nst->mins[i]) *(factor -1) *0.5;
  109.     nst->mins[i] -= t;          /* adapt the minimal */
  110.     nst->maxs[i] += t;          /* and   the maximal value */
  111.   }                             /* for all dimensions in range */
  112. }  /* nst_expand() */
  113. /*--------------------------------------------------------------------*/
  114. void nst_scale (NSTATS *nst, int idx, double off, double fac)
  115. {                               /* --- set (linear) scaling */
  116.   int i;                        /* loop variable */
  117.   assert(nst && (idx < nst->dim));  /* check the arguments */
  118.   if (idx < 0) { i = nst->dim; idx = 0; }
  119.   else         { i = idx +1; }  /* get index range to set */
  120.   while (--i >= idx) {          /* and traverse it */
  121.     nst->offs[i] = off;         /* set the offset */
  122.     nst->facs[i] = fac;         /* and the scaling factor */
  123.   }                             /* for all dimensions in range */
  124. }  /* nst_scale() */
  125. /*--------------------------------------------------------------------*/
  126. void nst_norm (NSTATS *nst, const double *vec, double *res)
  127. {                               /* --- normalize a data vector */
  128.   int    i;                     /* loop variable */
  129.   double *off, *fac;            /* to traverse the scaling parameters */
  130.   assert(nst && vec && res);    /* check the function arguments */
  131.   off = nst->offs +(i = nst->dim);
  132.   fac = nst->facs + i;          /* get the scaling parameters */
  133.   res += i; vec += i;           /* and the data vectors */
  134.   while (--i >= 0) *--res = *--fac * (*--vec - *--off);
  135. }  /* nst_norm() */             /* scale the vector */
  136. /*--------------------------------------------------------------------*/
  137. void nst_inorm (NSTATS *nst, const double *vec, double *res)
  138. {                               /* --- inverse normalize a vector */
  139.   int    i;                     /* loop variable */
  140.   double *off, *fac;            /* to traverse the scaling parameters */
  141.   assert(nst && vec && res);    /* check the function arguments */
  142.   off = nst->offs +(i = nst->dim);
  143.   fac = nst->facs + i;          /* get the scaling parameters */
  144.   res += i; vec += i;           /* and the data vectors */
  145.   while (--i >= 0) *--res = *--vec / *--fac + *--off;
  146. }  /* nst_inorm() */            /* scale the vector */
  147. /*--------------------------------------------------------------------*/
  148. void nst_center (NSTATS *nst, double *vec)
  149. {                               /* --- get center of data space */
  150.   int    i;                     /* loop variable */
  151.   double *min, *max;            /* to traverse the ranges */
  152.   assert(nst && vec);           /* check the function arguments */
  153.   min = nst->mins;              /* get the range variables, */
  154.   max = nst->maxs;              /* traverse the dimensions, */
  155.   for (i = nst->dim; --i >= 0;) /* and compute the center vector */
  156.     vec[i] = 0.5 *(max[i] +min[i]);
  157. }  /* nst_center() */
  158. /*--------------------------------------------------------------------*/
  159. void nst_spans (NSTATS *nst, double *vec)
  160. {                               /* --- get spans of dimensions */
  161.   int    i;                     /* loop variable */
  162.   double *min, *max;            /* to traverse the ranges */
  163.   assert(nst && vec);           /* check the function arguments */
  164.   min = nst->mins;            
  165.   max = nst->maxs;              /* get the range variables, */ 
  166.   for (i = nst->dim; --i >= 0;) /* traverse the dimensions, */ 
  167.     vec[i] = max[i] -min[i];    /* and compute the spans */
  168. }  /* nst_spans() */
  169. /*--------------------------------------------------------------------*/
  170. int nst_desc (NSTATS *nst, FILE *file, const char *indent, int maxlen)
  171. {                               /* --- describe norm. statistics */
  172.   int  i;                       /* loop variable */
  173.   int  pos, ind;                /* position in output line */
  174.   char buf[64];                 /* buffer for output */
  175.   for (i = nst->dim; --i >= 0;) /* check for non-identity scaling */
  176.     if ((nst->offs[i] != 0) || (nst->facs[i] != 1)) break;
  177.   if (i < 0) return 0;          /* if all identity scaling, abort */
  178.   fputs(indent,        file);   /* write the indentation and */
  179.   fputs("scales   = ", file);   /* start the scaling parameters */
  180.   for (ind = 0; indent[ind]; ind++);
  181.   pos = ind +9;                 /* compute the starting position */
  182.   for (i = 0; i < nst->dim; i++) {
  183.     pos += sprintf(buf, "[% g, %g]", nst->offs[i], nst->facs[i]);
  184.     if (i > 0) {                /* format the scaling parameters */
  185.       if (pos +3 <= maxlen) { fputs(", ", file);      pos += 2;   }
  186.       else { fprintf(file, ",n%s         ", indent); pos  = ind; }
  187.     }                           /* print separator and indentation */
  188.     fputs(buf, file);           /* print formatted offset and factor */
  189.   }
  190.   fputs(";n", file);           /* terminate the list */
  191.   return ferror(file) ? -1 : 0; /* return the write status */
  192. }  /* nst_desc() */
  193. /*--------------------------------------------------------------------*/
  194. #ifdef NST_PARSE
  195. static int _parse (SCAN *scan, int dim, double **buf)
  196. {                               /* --- parse normalization statistics */
  197.   int    k, n = 0;              /* loop variable, counter */
  198.   double *p;                    /* to access the statistics elements */
  199.   assert(scan);                 /* check the function arguments */
  200.   if ((sc_token(scan) != T_ID)  /* check whether 'scales' follows */
  201.   ||  (strcmp(sc_value(scan), "scales") != 0))
  202.     ERR_STR("scales");          /* if not, abort the function */
  203.   GET_TOK();                    /* consume 'scales' */
  204.   GET_CHR('=');                 /* consume '=' */
  205.   for (k = 0; (dim <= 0) || (k < dim); k++) {
  206.     if (k > 0) { GET_CHR(',');} /* if not first, consume ',' */
  207.     if (k >= n) {               /* if the statistics vector is full */
  208.       if (dim > 0) n  = dim;    /* compute the new vector size */
  209.       else         n += (n > BLKSIZE) ? n >> 1 : BLKSIZE;
  210.       p = (double*)realloc(*buf, (n+n) *sizeof(double));
  211.       if (!p) ERROR(E_NOMEM);   /* enlarge the buffer vector */
  212.       *buf = p;                 /* and set the new vector, */
  213.     }                           /* then note factor and offset */
  214.     p = *buf +k +k;             /* get the element to set */
  215.     GET_CHR('[');               /* consume '[' */
  216.     if (sc_token(scan) != T_NUM) ERROR(E_NUMEXP);
  217.     p[0] = strtod(sc_value(scan), NULL);
  218.     GET_TOK();                  /* consume the offset */
  219.     GET_CHR(',');               /* consume '[' */
  220.     if (sc_token(scan) != T_NUM) ERROR(E_NUMEXP);
  221.     p[1] = strtod(sc_value(scan), NULL);
  222.     GET_TOK();                  /* consume the factor */
  223.     GET_CHR(']');               /* consume '[' */
  224.     if ((dim <= 0) && (sc_token(scan) != ',')) {
  225.       k++; break; }             /* check for more scaling params. */
  226.   }
  227.   GET_CHR(';');                 /* consume ';' */
  228.   return k;                     /* return 'ok' */
  229. }  /* _parse() */
  230. /*--------------------------------------------------------------------*/
  231. NSTATS* nst_parse (SCAN *scan, int dim)
  232. {                               /* --- parse normalization statistics */
  233.   NSTATS *nst;                  /* created normalization statistics */
  234.   double *buf = NULL;           /* buffer for reading */
  235.   assert(scan);                 /* check the function arguments */
  236.   dim = _parse(scan,dim, &buf); /* parse normalization statistics */
  237.   if (dim < 0) { if (buf) free(buf); return NULL; }
  238.   nst = nst_create(dim);        /* create a statistics structure */
  239.   if (!nst)    { free(buf); return NULL; }
  240.   for (buf += dim +dim; --dim >= 0; ) {
  241.     nst->facs[dim] = *--buf;    /* copy the buffered values */
  242.     nst->offs[dim] = *--buf;    /* into the corresponding vectors */
  243.   }
  244.   free(buf);                    /* delete the read buffer */
  245.   return nst;                   /* return the created structure */
  246. }  /* nst_parse() */
  247. #endif