Sort all includes in src/gromacs
[alexxy/gromacs.git] / src / gromacs / gmxana / gmx_wham.cpp
index 7119bf4723612c99ded5174d062a7e81342e3a33..7a5b0a4333403b67af100b7fd78fb35d2466949f 100644 (file)
@@ -1,37 +1,38 @@
-/* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*- */
 /*
+ * This file is part of the GROMACS molecular simulation package.
  *
- *                This source code is part of
- *
- *                 G   R   O   M   A   C   S
- *
- *          GROningen MAchine for Chemical Simulations
- *
- *                        VERSION 3.2.0
- * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team,
- * check out http://www.gromacs.org for more information.
-
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
+ * Copyright (c) 2001-2004, The GROMACS development team.
+ * Copyright (c) 2013,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
  * of the License, or (at your option) any later version.
  *
- * If you want to redistribute modifications, please consider that
- * scientific software is very special. Version control is crucial -
- * bugs must be traceable. We will be happy to consider code for
- * inclusion in the official distribution, but derived work must not
- * be called official GROMACS. Details are found in the README & COPYING
- * files - if they are missing, get the official version at www.gromacs.org.
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
  *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the papers on the package - you can find them in the top README file.
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  *
- * For more info, check our website at http://www.gromacs.org
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
  *
- * And Hey:
- * Green Red Orange Magenta Azure Cyan Skyblue
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
  */
 
 /*! \internal \file
  *  \author Jochen Hub <jhub@gwdg.de>
  */
 
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
+#include "gmxpre.h"
+
+#include "config.h"
+
+#include <ctype.h>
 #include <stdio.h>
 #include <stdlib.h>
+#include <string.h>
+
 #include <sstream>
 
 #include "gromacs/commandline/pargs.h"
-#include "typedefs.h"
-#include "smalloc.h"
-#include "vec.h"
-#include "copyrite.h"
 #include "gromacs/fileio/tpxio.h"
-#include "names.h"
-#include "gmx_random.h"
-#include "gmx_ana.h"
-#include "macros.h"
-
-#include "string2.h"
-#include "xvgr.h"
+#include "gromacs/fileio/xvgr.h"
+#include "gromacs/gmxana/gmx_ana.h"
+#include "gromacs/legacyheaders/copyrite.h"
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/legacyheaders/names.h"
+#include "gromacs/legacyheaders/typedefs.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/random/random.h"
+#include "gromacs/utility/cstringutil.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/smalloc.h"
 
 //! longest file names allowed in input files
 #define WHAM_MAXFILELEN 2048
 
+/*! \brief
+ * x-axis legend for output files
+ */
+static const char *xlabel = "\\xx\\f{} (nm)";
+
 /*! \brief
  * enum for energy units
  */
@@ -110,14 +120,14 @@ typedef struct
      * \name Using umbrella pull code since gromacs 4.x
      */
     /*!\{*/
-    int   npullcrds;     //!< nr of pull coordinates in tpr file
-    int   pull_geometry; //!< such as distance, direction
-    ivec  pull_dim;      //!< pull dimension with geometry distance
-    int   pull_ndim;     //!< nr of pull_dim != 0
-    gmx_bool bPrintRef;  //!< Coordinates of reference groups written to pullx.xvg ?
-    real *k;             //!< force constants in tpr file
-    real *init_dist;     //!< reference displacements
-    real *umbInitDist;   //!< reference displacement in umbrella direction
+    int      npullcrds;     //!< nr of pull coordinates in tpr file
+    int      pull_geometry; //!< such as distance, direction
+    ivec     pull_dim;      //!< pull dimension with geometry distance
+    int      pull_ndim;     //!< nr of pull_dim != 0
+    gmx_bool bPrintRef;     //!< Coordinates of reference groups written to pullx.xvg ?
+    real    *k;             //!< force constants in tpr file
+    real    *init_dist;     //!< reference displacements
+    real    *umbInitDist;   //!< reference displacement in umbrella direction
     /*!\}*/
     /*!
      * \name Using PDO files common until gromacs 3.x
@@ -1203,7 +1213,7 @@ void calc_cumulatives(t_UmbrellaWindow *window, int nWindows,
         snew(fn, strlen(fnhist)+10);
         snew(buf, strlen(fnhist)+10);
         sprintf(fn, "%s_cumul.xvg", strncpy(buf, fnhist, strlen(fnhist)-4));
-        fp = xvgropen(fn, "CDFs of umbrella windows", "z", "CDF", opt->oenv);
+        fp = xvgropen(fn, "CDFs of umbrella windows", xlabel, "CDF", opt->oenv);
     }
 
     nbin = opt->bins;
@@ -1244,7 +1254,7 @@ void calc_cumulatives(t_UmbrellaWindow *window, int nWindows,
             fprintf(fp, "\n");
         }
         printf("Wrote cumulative distribution functions to %s\n", fn);
-        ffclose(fp);
+        gmx_ffclose(fp);
         sfree(fn);
         sfree(buf);
     }
@@ -1449,11 +1459,11 @@ void print_histograms(const char *fnhist, t_UmbrellaWindow * window, int nWindow
     }
     else
     {
-        fn = strdup(fnhist);
+        fn = gmx_strdup(fnhist);
         strcpy(title, "Umbrella histograms");
     }
 
-    fp   = xvgropen(fn, title, "z", "count", opt->oenv);
+    fp   = xvgropen(fn, title, xlabel, "count", opt->oenv);
     bins = opt->bins;
 
     /* Write histograms */
@@ -1470,7 +1480,7 @@ void print_histograms(const char *fnhist, t_UmbrellaWindow * window, int nWindow
         fprintf(fp, "\n");
     }
 
-    ffclose(fp);
+    gmx_ffclose(fp);
     printf("Wrote %s\n", fn);
     if (bs_index >= 0)
     {
@@ -1626,7 +1636,7 @@ void do_bootstrapping(const char *fnres, const char* fnprof, const char *fnhist,
     }
 
     /* do bootstrapping */
-    fp = xvgropen(fnprof, "Boot strap profiles", "z", ylabel, opt->oenv);
+    fp = xvgropen(fnprof, "Boot strap profiles", xlabel, ylabel, opt->oenv);
     for (ib = 0; ib < opt->nBootStrap; ib++)
     {
         printf("  *******************************************\n"
@@ -1712,13 +1722,16 @@ void do_bootstrapping(const char *fnres, const char* fnprof, const char *fnhist,
             bsProfiles_av2[i] += tmp*tmp;
             fprintf(fp, "%e\t%e\n", (i+0.5)*opt->dz+opt->min, tmp);
         }
-        fprintf(fp, "&\n");
+        fprintf(fp, "%s\n", output_env_get_print_xvgr_codes(opt->oenv) ? "&" : "");
     }
-    ffclose(fp);
+    gmx_ffclose(fp);
 
     /* write average and stddev */
-    fp = xvgropen(fnres, "Average and stddev from bootstrapping", "z", ylabel, opt->oenv);
-    fprintf(fp, "@TYPE xydy\n");
+    fp = xvgropen(fnres, "Average and stddev from bootstrapping", xlabel, ylabel, opt->oenv);
+    if (output_env_get_print_xvgr_codes(opt->oenv))
+    {
+        fprintf(fp, "@TYPE xydy\n");
+    }
     for (i = 0; i < opt->bins; i++)
     {
         bsProfiles_av [i] /= opt->nBootStrap;
@@ -1727,7 +1740,7 @@ void do_bootstrapping(const char *fnres, const char* fnprof, const char *fnhist,
         stddev             = (tmp >= 0.) ? sqrt(tmp) : 0.; /* Catch rouding errors */
         fprintf(fp, "%e\t%e\t%e\n", (i+0.5)*opt->dz+opt->min, bsProfiles_av [i], stddev);
     }
-    ffclose(fp);
+    gmx_ffclose(fp);
     printf("Wrote boot strap result to %s\n", fnres);
 }
 
@@ -1763,7 +1776,7 @@ void read_wham_in(const char *fn, char ***filenamesRet, int *nfilesRet,
     int    nread, sizenow, i, block = 1;
     FILE  *fp;
 
-    fp      = ffopen(fn, "r");
+    fp      = gmx_ffopen(fn, "r");
     nread   = 0;
     sizenow = 0;
     while (fgets(tmp, sizeof(tmp), fp) != NULL)
@@ -1860,7 +1873,7 @@ FILE *open_pdo_pipe(const char *fn, t_UmbrellaOptions *opt, gmx_bool *bPipeOpen)
     }
     else
     {
-        pipe       = ffopen(fn, "r");
+        pipe       = gmx_ffopen(fn, "r");
         *bPipeOpen = FALSE;
     }
 
@@ -1873,7 +1886,7 @@ void pdo_close_file(FILE *fp)
 #ifdef HAVE_PIPES
     pclose(fp);
 #else
-    ffclose(fp);
+    gmx_ffclose(fp);
 #endif
 }
 
@@ -1926,7 +1939,7 @@ void read_pdo_files(char **fn, int nfiles, t_UmbrellaHeader* header,
             }
             else
             {
-                ffclose(file);
+                gmx_ffclose(file);
             }
         }
         printf("\n");
@@ -1964,7 +1977,7 @@ void read_pdo_files(char **fn, int nfiles, t_UmbrellaHeader* header,
         }
         else
         {
-            ffclose(file);
+            gmx_ffclose(file);
         }
     }
     printf("\n");
@@ -2679,7 +2692,7 @@ void calcIntegratedAutocorrelationTimes(t_UmbrellaWindow *window, int nwins,
                 {
                     fprintf(fpcorr, "%g  %g\n", k*dt, corr[k]);
                 }
-                fprintf(fpcorr, "&\n");
+                fprintf(fpcorr, "%s\n", output_env_get_print_xvgr_codes(opt->oenv) ? "&" : "");
             }
 
             /* esimate integrated correlation time, fitting is too unstable */
@@ -2702,21 +2715,24 @@ void calcIntegratedAutocorrelationTimes(t_UmbrellaWindow *window, int nwins,
     printf(" done\n");
     if (fpcorr)
     {
-        ffclose(fpcorr);
+        gmx_ffclose(fpcorr);
     }
 
     /* plot IACT along reaction coordinate */
-    fp = xvgropen(fn, "Integrated autocorrelation times", "z", "IACT [ps]", opt->oenv);
-    fprintf(fp, "@    s0 symbol 1\n@    s0 symbol size 0.5\n@    s0 line linestyle 0\n");
-    fprintf(fp, "#  WIN   tau(gr1)  tau(gr2) ...\n");
-    for (i = 0; i < nwins; i++)
+    fp = xvgropen(fn, "Integrated autocorrelation times", xlabel, "IACT [ps]", opt->oenv);
+    if (output_env_get_print_xvgr_codes(opt->oenv))
     {
-        fprintf(fp, "# %3d   ", i);
-        for (ig = 0; ig < window[i].nPull; ig++)
+        fprintf(fp, "@    s0 symbol 1\n@    s0 symbol size 0.5\n@    s0 line linestyle 0\n");
+        fprintf(fp, "#  WIN   tau(gr1)  tau(gr2) ...\n");
+        for (i = 0; i < nwins; i++)
         {
-            fprintf(fp, " %11g", window[i].tau[ig]);
+            fprintf(fp, "# %3d   ", i);
+            for (ig = 0; ig < window[i].nPull; ig++)
+            {
+                fprintf(fp, " %11g", window[i].tau[ig]);
+            }
+            fprintf(fp, "\n");
         }
-        fprintf(fp, "\n");
     }
     for (i = 0; i < nwins; i++)
     {
@@ -2731,8 +2747,12 @@ void calcIntegratedAutocorrelationTimes(t_UmbrellaWindow *window, int nwins,
                opt->sigSmoothIact);
         /* smooth IACT along reaction coordinate and overwrite g=1+2tau */
         smoothIact(window, nwins, opt);
-        fprintf(fp, "&\n@    s1 symbol 1\n@    s1 symbol size 0.5\n@    s1 line linestyle 0\n");
-        fprintf(fp, "@    s1 symbol color 2\n");
+        fprintf(fp, "%s\n", output_env_get_print_xvgr_codes(opt->oenv) ? "&" : "");
+        if (output_env_get_print_xvgr_codes(opt->oenv))
+        {
+            fprintf(fp, "@    s1 symbol 1\n@    s1 symbol size 0.5\n@    s1 line linestyle 0\n");
+            fprintf(fp, "@    s1 symbol color 2\n");
+        }
         for (i = 0; i < nwins; i++)
         {
             for (ig = 0; ig < window[i].nPull; ig++)
@@ -2741,7 +2761,7 @@ void calcIntegratedAutocorrelationTimes(t_UmbrellaWindow *window, int nwins,
             }
         }
     }
-    ffclose(fp);
+    gmx_ffclose(fp);
     printf("Wrote %s\n", fn);
 }
 
@@ -2939,7 +2959,7 @@ void guessPotByIntegration(t_UmbrellaWindow *window, int nWindows, t_UmbrellaOpt
                     nHist++;
                     fAv += window[i].forceAv[ig];
                 }
-                /* at the same time, rememer closest histogram */
+                /* at the same time, remember closest histogram */
                 if (dist < distmin)
                 {
                     winmin   = i;
@@ -2975,12 +2995,12 @@ void guessPotByIntegration(t_UmbrellaWindow *window, int nWindows, t_UmbrellaOpt
     }
     if (opt->verbose)
     {
-        fp = xvgropen("pmfintegrated.xvg", "PMF from force integration", "z", "PMF [kJ/mol]", opt->oenv);
+        fp = xvgropen("pmfintegrated.xvg", "PMF from force integration", xlabel, "PMF (kJ/mol)", opt->oenv);
         for (j = 0; j < opt->bins; ++j)
         {
             fprintf(fp, "%g  %g\n", (j+0.5)*dz+opt->min, pot[j]);
         }
-        ffclose(fp);
+        gmx_ffclose(fp);
         printf("verbose mode: wrote %s with PMF from interated forces\n", "pmfintegrated.xvg");
     }
 
@@ -3034,7 +3054,7 @@ void readPullGroupSelection(t_UmbrellaOptions *opt, char **fnTpr, int nTpr)
     char  fmt[1024], fmtign[1024];
     int   block = 1, sizenow;
 
-    fp            = ffopen(opt->fnGroupsel, "r");
+    fp            = gmx_ffopen(opt->fnGroupsel, "r");
     opt->groupsel = NULL;
 
     snew(tmpbuf, len);
@@ -3099,7 +3119,7 @@ void readPullGroupSelection(t_UmbrellaOptions *opt, char **fnTpr, int nTpr)
 //! Boolean XOR
 #define WHAMBOOLXOR(a, b) ( ((!(a)) && (b)) || ((a) && (!(b))))
 
-/*! Number of elements in fnm (used for command line parsing) */
+//! Number of elements in fnm (used for command line parsing)
 #define NFILE asize(fnm)
 
 //! The main g_wham routine
@@ -3365,7 +3385,7 @@ int gmx_wham(int argc, char *argv[])
     opt.stepchange            = 100;
     opt.stepUpdateContrib     = 100;
 
-    if (!parse_common_args(&argc, argv, PCA_BE_NICE,
+    if (!parse_common_args(&argc, argv, 0,
                            NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &opt.oenv))
     {
         return 0;
@@ -3486,7 +3506,7 @@ int gmx_wham(int argc, char *argv[])
 
     /* write histograms */
     histout = xvgropen(opt2fn("-hist", NFILE, fnm), "Umbrella histograms",
-                       "z", "count", opt.oenv);
+                       xlabel, "count", opt.oenv);
     for (l = 0; l < opt.bins; ++l)
     {
         fprintf(histout, "%e\t", (double)(l+0.5)/opt.bins*(opt.max-opt.min)+opt.min);
@@ -3499,7 +3519,7 @@ int gmx_wham(int argc, char *argv[])
         }
         fprintf(histout, "\n");
     }
-    ffclose(histout);
+    gmx_ffclose(histout);
     printf("Wrote %s\n", opt2fn("-hist", NFILE, fnm));
     if (opt.bHistOnly)
     {
@@ -3595,12 +3615,12 @@ int gmx_wham(int argc, char *argv[])
     }
 
     /* write profile or density of states */
-    profout = xvgropen(opt2fn("-o", NFILE, fnm), title, "z", ylabel, opt.oenv);
+    profout = xvgropen(opt2fn("-o", NFILE, fnm), title, xlabel, ylabel, opt.oenv);
     for (i = 0; i < opt.bins; ++i)
     {
         fprintf(profout, "%e\t%e\n", (double)(i+0.5)/opt.bins*(opt.max-opt.min)+opt.min, profile[i]);
     }
-    ffclose(profout);
+    gmx_ffclose(profout);
     printf("Wrote %s\n", opt2fn("-o", NFILE, fnm));
 
     /* Bootstrap Method */