src/gromacs/legacyheaders/gmx_statistics.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2008, The GROMACS development team.
   6  * Copyright (c) 2010, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37
  38 #ifndef _GMX_STATS_H
  39 #define _GMX_STATS_H
  40
  41 #ifdef __cplusplus
  42 extern "C" {
  43 #endif
  44
  45 #include "typedefs.h"
  46
  47 typedef struct gmx_stats *gmx_stats_t;
  48
  49 /* Error codes returned by the routines */
  50 enum {
  51     estatsOK, estatsNO_POINTS, estatsNO_MEMORY, estatsERROR,
  52     estatsINVALID_INPUT, estatsNOT_IMPLEMENTED, estatsNR
  53 };
  54
  55 enum {
  56     elsqWEIGHT_NONE, elsqWEIGHT_X, elsqWEIGHT_Y,
  57     elsqWEIGHT_XY, elsqWEIGHT_NR
  58 };
  59
  60 enum {
  61     ehistoX, ehistoY, ehistoNR
  62 };
  63
  64 gmx_stats_t gmx_stats_init();
  65
  66 int gmx_stats_done(gmx_stats_t stats);
  67
  68 /* Remove outliers from a straight line, where level in units of
  69    sigma. Level needs to be larger than one obviously. */
  70 int gmx_stats_remove_outliers(gmx_stats_t stats, double level);
  71
  72
  73
  74 int gmx_stats_add_point(gmx_stats_t stats, double x, double y,
  75                         double dx, double dy);
  76
  77 /* The arrays dx and dy may be NULL if no uncertainties are available,
  78    in that case zero uncertainties will be assumed. */
  79 int gmx_stats_add_points(gmx_stats_t stats, int n, real *x, real *y,
  80                          real *dx, real *dy);
  81
  82 /* Return the data points one by one. Return estatsOK while there are
  83    more points, and returns estatsNOPOINTS when the last point has
  84    been returned. Should be used in a while loop. Variables for either
  85    pointer may be NULL, in which case the routine can be used as an
  86    expensive point counter.
  87    If level > 0 then the outliers outside level*sigma are reported
  88    only. */
  89 int gmx_stats_get_point(gmx_stats_t stats, real *x, real *y,
  90                         real *dx, real *dy, real level);
  91
  92 /* Fit the data to y = ax + b, possibly weighted, if uncertainties
  93    have been input. Returns slope in *a and intercept in b, *return
  94    sigmas in *da and *db respectively. Returns normalized *quality of
  95    fit in *chi2 and correlation of fit with data in Rfit. chi2, Rfit,
  96    da and db may be NULL. */
  97 int gmx_stats_get_ab(gmx_stats_t stats, int weight,
  98                      real *a, real *b,
  99                      real *da, real *db, real *chi2, real *Rfit);
 100
 101 /* Fit the data to y = ax, possibly weighted, if uncertainties have
 102    been input. Returns slope in *a, sigma in a in *da, and normalized
 103    quality of fit in *chi2 and correlation of fit with data in
 104    Rfit. chi2, Rfit and da may be NULL. */
 105 int gmx_stats_get_a(gmx_stats_t stats, int weight,
 106                     real *a, real *da, real *chi2, real *Rfit);
 107
 108 /* Return the correlation coefficient between the data (x and y) as
 109    input to the structure. */
 110 int gmx_stats_get_corr_coeff(gmx_stats_t stats, real *R);
 111
 112 /* Returns the root mean square deviation between x and y values. */
 113 int gmx_stats_get_rmsd(gmx_stats_t gstats, real *rmsd);
 114
 115 int gmx_stats_get_npoints(gmx_stats_t stats, int *N);
 116
 117 int gmx_stats_get_average(gmx_stats_t stats, real *aver);
 118
 119 int gmx_stats_get_sigma(gmx_stats_t stats, real *sigma);
 120
 121 int gmx_stats_get_error(gmx_stats_t stats, real *error);
 122
 123 /* Get all three of the above. Pointers may be null, in which case no
 124    assignment will be done. */
 125 int gmx_stats_get_ase(gmx_stats_t gstats, real *aver, real *sigma, real *error);
 126
 127 /* Dump the x, y, dx, dy data to a text file */
 128 int gmx_stats_dump_xy(gmx_stats_t gstats, FILE *fp);
 129
 130 /* Make a histogram of the data present. Uses either bindwith to
 131    determine the number of bins, or nbins to determine the binwidth,
 132    therefore one of these should be zero, but not the other. If *nbins = 0
 133    the number of bins will be returned in this variable. ehisto should be one of
 134    ehistoX or ehistoY. If
 135    normalized not equal to zero, the integral of the histogram will be
 136    normalized to one. The output is in two arrays, *x and *y, to which
 137    you should pass a pointer. Memory for the arrays will be allocated
 138    as needed. Function returns one of the estats codes. */
 139 int gmx_stats_make_histogram(gmx_stats_t gstats, real binwidth, int *nbins,
 140                              int ehisto,
 141                              int normalized, real **x, real **y);
 142
 143 /* Return message belonging to error code */
 144 const char *gmx_stats_message(int estats);
 145
 146 /****************************************************
 147  * Some statistics utilities for convenience: useful when a complete data
 148  * set is available already from another source, e.g. an xvg file.
 149  ****************************************************/
 150 int lsq_y_ax(int n, real x[], real y[], real *a);
 151 /* Fit a straight line y=ax thru the n data points x, y, return the
 152    slope in *a. Return value can be estatsOK, or something else. */
 153
 154 int lsq_y_ax_b(int n, real x[], real y[], real *a, real *b, real *r,
 155                real *chi2);
 156 /* Fit a straight line y=ax+b thru the n data points x,y.
 157  * Returns the "fit quality" sigma = sqrt(chi^2/(n-2)).
 158  * The correlation coefficient is returned in r.
 159  */
 160
 161 int lsq_y_ax_b_xdouble(int n, double x[], real y[],
 162                        real *a, real *b, real *r, real *chi2);
 163 /* As lsq_y_ax_b, but with x in double precision.
 164  */
 165
 166 int lsq_y_ax_b_error(int n, real x[], real y[], real dy[],
 167                      real *a, real *b, real *da, real *db,
 168                      real *r, real *chi2);
 169 /* Fit a straight line y=ax+b thru the n data points x,y, with sigma dy
 170  * Returns the "fit quality" sigma = sqrt(chi^2/(n-2)).
 171  * The correlation coefficient is returned in r.
 172  */
 173
 174 #ifdef __cplusplus
 175 }
 176 #endif
 177
 178 #endif