include/gmx_statistics.h

   1 #ifndef _GMX_STATS_H
   2 #define _GMX_STATS_H
   3
   4 typedef struct gmx_stats *gmx_stats_t;
   5
   6 /* Error codes returned by the routines */
   7 enum { estatsOK, estatsNO_POINTS, estatsNO_MEMORY, estatsERROR,
   8        estatsINVALID_INPUT, estatsNOT_IMPLEMENTED, estatsNR };
   9
  10 enum { elsqWEIGHT_NONE, elsqWEIGHT_X, elsqWEIGHT_Y,
  11        elsqWEIGHT_XY, elsqWEIGHT_NR  };
  12
  13 extern gmx_stats_t gmx_stats_init();
  14
  15 extern int gmx_stats_done(gmx_stats_t stats);
  16
  17 /* Remove outliers from a straight line, where level in units of
  18    sigma. Level needs to be larger than one obviously. */
  19 extern int gmx_stats_remove_outliers(gmx_stats_t stats,double level);
  20
  21 extern int gmx_stats_add_point(gmx_stats_t stats,double x,double y,
  22                                double dx,double dy);
  23
  24 /* The arrays dx and dy may be NULL if no uncertainties are available,
  25    in that case zero uncertainties will be assumed. */
  26 extern int gmx_stats_add_points(gmx_stats_t stats,int n,real *x,real *y,
  27                                 real *dx,real *dy);
  28
  29 /* Return the data points one by one. Return estatsOK while there are
  30    more points, and returns estatsNOPOINTS when the last point has
  31    been returned. Should be used in a while loop. Variables for either
  32    pointer may be NULL, in which case the routine can be used as an
  33    expensive point counter. */
  34 extern int gmx_stats_get_point(gmx_stats_t stats,real *x,real *y,
  35                                real *dx,real *dy);
  36
  37 /* Fit the data to y = ax + b, possibly weighted, if uncertainties
  38    have been input. Returns slope in *a and intercept in b, *return
  39    sigmas in *da and *db respectively. Returns normalized *quality of
  40    fit in *chi2 and correlation of fit with data in Rfit. chi2, Rfit,
  41    da and db may be NULL. */
  42 extern int gmx_stats_get_ab(gmx_stats_t stats,int weight,
  43                             real *a,real *b,
  44                             real *da,real *db,real *chi2,real *Rfit);
  45
  46 /* Fit the data to y = ax, possibly weighted, if uncertainties have
  47    been input. Returns slope in *a, sigma in a in *da, and normalized
  48    quality of fit in *chi2 and correlation of fit with data in
  49    Rfit. chi2, Rfit and da may be NULL. */
  50 extern int gmx_stats_get_a(gmx_stats_t stats,int weight,
  51                            real *a,real *da,real *chi2,real *Rfit);
  52
  53 /* Return the correlation coefficient between the data (x and y) as
  54    input to the structure. */
  55 extern int gmx_stats_get_corr_coeff(gmx_stats_t stats,real *R);
  56
  57 /* Returns the root mean square deviation between x and y values. */
  58 extern int gmx_stats_get_rmsd(gmx_stats_t gstats,real *rmsd);
  59
  60 extern int gmx_stats_get_npoints(gmx_stats_t stats,int *N);
  61
  62 extern int gmx_stats_get_average(gmx_stats_t stats,real *aver);
  63
  64 extern int gmx_stats_get_sigma(gmx_stats_t stats,real *sigma);
  65
  66 extern int gmx_stats_get_error(gmx_stats_t stats,real *error);
  67
  68 /* Get all three of the above. Pointers may be null, in which case no
  69    assignment will be done. */
  70 extern int gmx_stats_get_ase(gmx_stats_t gstats,real *aver,real *sigma,real *error);
  71
  72 /* Dump the x, y, dx, dy data to a text file */
  73 extern int gmx_stats_dump_xy(gmx_stats_t gstats,FILE *fp);
  74
  75 /* Make a histogram of the data present. Uses either bindwith to
  76    determine the number of bins, or nbins to determine the binwidth,
  77    therefore one of these should be zero, but not the other. If
  78    normalized not equal to zero, the integral of the histogram will be
  79    normalized to one. The output is in two arrays, *x and *y, to which
  80    you should pass a pointer. Memory for the arrays will be allocated
  81    as needed. Function returns one of the estats codes. */
  82 extern int gmx_stats_make_histogram(gmx_stats_t gstats,real binwidth,int nbins,
  83                                     int normalized,real **x,real **y);
  84
  85 /* Return message belonging to error code */
  86 extern const char *gmx_stats_message(int estats);
  87
  88 /****************************************************
  89  * Some statistics utilities for convenience: useful when a complete data
  90  * set is available already from another source, e.g. an xvg file.
  91  ****************************************************/
  92 extern int lsq_y_ax(int n, real x[], real y[], real *a);
  93 /* Fit a straight line y=ax thru the n data points x, y, return the
  94    slope in *a. Return value can be estatsOK, or something else. */
  95
  96 extern int lsq_y_ax_b(int n, real x[], real y[], real *a, real *b,real *r,
  97                       real *chi2);
  98 /* Fit a straight line y=ax+b thru the n data points x,y.
  99  * Returns the "fit quality" sigma = sqrt(chi^2/(n-2)).
 100  * The correlation coefficient is returned in r.
 101  */
 102
 103 extern int lsq_y_ax_b_xdouble(int n, double x[], real y[],
 104                                real *a, real *b,real *r,real *chi2);
 105 /* As lsq_y_ax_b, but with x in double precision.
 106  */
 107
 108 extern int lsq_y_ax_b_error(int n, real x[], real y[], real dy[],
 109                              real *a, real *b, real *da, real *db,
 110                              real *r,real *chi2);
 111 /* Fit a straight line y=ax+b thru the n data points x,y, with sigma dy
 112  * Returns the "fit quality" sigma = sqrt(chi^2/(n-2)).
 113  * The correlation coefficient is returned in r.
 114  */
 115
 116 #endif